home *** CD-ROM | disk | FTP | other *** search
/ Aminet 35 / Aminet 35 (2000)(Schatztruhe)[!][Feb 2000].iso / Aminet / gfx / misc / gnuplot-src.lha / gnuplot-3.7.1src / gnuplot-3.7.1.lha / gnuplot-3.7.1 / hidden3d.c < prev    next >
Encoding:
C/C++ Source or Header  |  1999-08-19  |  84.8 KB  |  2,796 lines
  1. #ifndef lintound repSid = "$Id: hidden3d.c,v 1.13.2.1 1999/08/19 14:39:28 lhecking Exp $";
  2. #endif
  3.  
  4. /* GNUPLOT - hidden3d.c */
  5.  
  6. /*[
  7.  * Copyright 1986 - 1993, 1998   Thomas Williams, Colin Kelley
  8.  *
  9.  * Permission to use, copy, and distribute this software and its
  10.  * documentation for any purpose with or without fee is hereby granted,
  11.  * provided that the above copyright notice appear in all copies and
  12.  * that both that copyright notice and this permission notice appear
  13.  * in supporting documentation.
  14.  *
  15.  * Permission to modify the software is granted, but not the right to
  16.  * distribute the complete modified source code.  Modifications are to
  17.  * be distributed as patches to the released version.  Permission to
  18.  * distribute binaries produced by compiling modified sources is granted,
  19.  * provided you
  20.  *   1. distribute the corresponding source modifications from the
  21.  *    released version in the form of a patch file along with the binaries,
  22.  *   2. add special version identification to distinguish your version
  23.  *    in addition to the base release version number,
  24.  *   3. provide your name and address as the primary contact for the
  25.  *    support of your modified version, and
  26.  *   4. retain our contact information in regard to use of the base
  27.  *    software.
  28.  * Permission to distribute the released version of the source code along
  29.  * with corresponding source modifications in the form of a patch file is
  30.  * granted with same provisions 2 through 4 for binary distributions.
  31.  *
  32.  * This software is provided "as is" without express or implied warranty
  33.  * to the extent permitted by applicable law.
  34. ]*/
  35.  
  36.  
  37. /*
  38.  * 19 September 1992  Lawrence Crowl  (crowl@cs.orst.edu)
  39.  * Added user-specified bases for log scaling.
  40.  *
  41.  * 'someone'  contributed a complete rewrite of the hidden line
  42.  * stuff, for about beta 173 or so, called 'newhide.tgz'
  43.  *
  44.  * 1995-1997 Hans-Bernhard Br"oker (Broeker@physik.rwth-aachen.de)
  45.  *   Speedup, cleanup and several modification to integrate
  46.  *   'newhide' with 3.6 betas ~188 (start of work) up to 273.
  47.  *   As of beta 290, this is 'officially' in gnuplot.
  48.  *
  49.  */
  50.  
  51. #include "plot.h"
  52. #include "setshow.h"
  53.  
  54. /* TODO (HBB's notes, just in case you're interested):
  55.  * + fix all the problems I annotated by a 'FIXME' comment
  56.  * + Find out which value EPSILON should have, and why
  57.  * + Ask gnuplot-beta for a suitable way of concentrating
  58.  *   all those 'VERYSMALL', 'EPSILON', and other numbers
  59.  *   into one, consistent scheme, possibly based on
  60.  *   <float.h>. E.g., I'd say that for most applications,
  61.  *   the best 'epsilon' is either DBL_EPS or its square root.
  62.  * + redo all the profiling, esp. to find out if TEST_GRIDCHECK = 1
  63.  *   actually gains some speed, and if it is worth the memory
  64.  *   spent to store the bit masks
  65.  *   -> seems not improve speed at all, at least for my standard
  66.  *   test case (mandd.gpl) it even slows it down by ~ 10%
  67.  * + Evaluate the speed/memory comparison for storing vertex
  68.  *   indices instead of (double) floating point constants
  69.  *   to store {x|y}{min|max}
  70.  * + remove those of my comments that are now pointless
  71.  * + indent all that code
  72.  * + try to really understand that 'hl_buffer' stuff...
  73.  * + get rid of hardcoded things like sizeof(short) == 4,
  74.  *   those 0x7fff terms and all that.
  75.  * + restructure the hl_buffer storing method to make it more efficient
  76.  * + Try to cut the speed decrease of this code rel. to the old hidden
  77.  *   hidden line removal. (For 'mandd.gpl', it costs an overall
  78.  *   factor of 9 compared with my older versions of gnuplot!)
  79.  */
  80.  
  81. /*************************/
  82. /* Configuration section */
  83. /*************************/
  84.  
  85. /* HBB: if this module is compiled with TEST_GRIDCHECK = 1 defined,
  86.  * it will store the information about {x|y}{min|max} in an
  87.  * other (additional!) form: a bit mask, with each bit representing
  88.  * one horizontal or vertical strip of the screen. The bits
  89.  * for  strips a polygon spans are set to one. This allows to
  90.  * test for xy_overlap simply by comparing bit patterns.
  91.  */
  92. #ifndef TEST_GRIDCHECK
  93. #define TEST_GRIDCHECK 0
  94. #endif
  95.  
  96. /* HBB 961124; If you don't want the color-distinction between the
  97.  * 'top' and 'bottom' sides of the surface, like I do, then just compile
  98.  * with -DBACKSIDE_LINETYPE_OFFSET = 0. */
  99. #ifndef BACKSIDE_LINETYPE_OFFSET
  100. #define BACKSIDE_LINETYPE_OFFSET 1
  101. #endif
  102.  
  103. /* HBB 961127: defining FOUR_TRIANGLES = 1 will separate each quadrangle
  104.  * of data points into *four* triangles, by introducing an additional
  105.  * point at the mean of the four corners. Status: experimental 
  106.  */
  107. #ifndef FOUR_TRIANGLES
  108. #define FOUR_TRIANGLES 0
  109. #endif
  110.  
  111. /* HBB 970322: I removed the SENTINEL flag and its implementation
  112.  * completely.  Its speed improvement wasn't that great, and it never
  113.  * fully worked anyway, so that shouldn't make much of a difference at
  114.  * all. */
  115.  
  116. /* HBB 961212: this #define lets you choose if the diagonals that
  117.  * divide each original quadrangle in two triangles will be drawn
  118.  * visible or not: do draw them, define it to be 7L, otherwise let be
  119.  * 3L (for the FOUR_TRIANGLES method, the values are 7L and 1L) */
  120. /* drd : default to 3, for back-compatibility with 3.5 */
  121. #ifndef TRIANGLE_LINESDRAWN_PATTERN
  122. #define TRIANGLE_LINESDRAWN_PATTERN 3L
  123. #endif
  124.  
  125. /* HBB 970131: with HANDLE_UNDEFINED_POINTS = 1, let's try to handle
  126.  * UNDEFINED data points in the input we're given by the rest of
  127.  * gnuplot. We do that by marking these points by giving them z = -2
  128.  * (which normally never happens), and then refusing to build any
  129.  * polygon if one of its vertices has this mark. Thus, there's now a
  130.  * hole in the generated surface. */
  131. /* HBB 970322: some of this code is now hardwired (in
  132.  * PREPARE_POLYGON), so HANDLE_UNDEFINED_POINTS = 0 might have stopped
  133.  * to be a useful setting. I doubt anyone will miss it, anyway. */
  134. /* drd : 2 means undefined only. 1 means outrange treated as undefined */
  135. #ifndef HANDLE_UNDEFINED_POINTS
  136. #define HANDLE_UNDEFINED_POINTS 1
  137. #endif
  138. /* Symbolic value for 'do not handle Undefined Points specially' */
  139. #define UNHANDLED (UNDEFINED+1)
  140.  
  141. /* HBB 970322: If both subtriangles of a quad were cancelled, try if
  142.  * using the other diagonal is better. This only makes a difference if
  143.  * exactly one vertex of the quad is unusable, and that one is on the
  144.  * first tried diagonal. In such a case, the border of the hole in the
  145.  * surface will be less rough than with the previous method. To
  146.  * disable this, #define SHOW_ALTERNATIVE_DIAGONAL as 0 */
  147. #ifndef SHOW_ALTERNATIVE_DIAGONAL
  148. #define SHOW_ALTERNATIVE_DIAGONAL 1
  149. #endif
  150.  
  151. /* HBB 9790522: If the two triangles in a quad are both drawn, and
  152.  * they show different sides to the user (the quad is 'bent over'),
  153.  * then it often looks more sensible to draw the diagonal visibly to
  154.  * avoid other parts of the scene being obscured by what the user
  155.  * can't see. To disable, #define HANDLE_BENTOVER_QUADRANGLES 0 */
  156. #ifndef HANDLE_BENTOVER_QUADRANGLES
  157. #define HANDLE_BENTOVER_QUADRANGLES 1
  158. #endif
  159.  
  160. /* HBB 970618: The code of split_polygon_by_plane used to make a
  161.  * separate decision about what should happen to points that are 'in'
  162.  * the splitting plane (within EPSILON accuracy, i.e.), based on the
  163.  * orientation of the splitting plane. I had disabled that code for I
  164.  * assumed it might be the cause of some of the buggyness. I'm not yet
  165.  * fully convinced on wether that assumption holds or not, so it's now
  166.  * choosable.  OLD_SPLIT_PLANE == 1 will enable it. Some older comments
  167.  * from the source:*/
  168. /* HBB 970325: re-inserted this from older versions. Finally, someone
  169.  * came up with a test case where it is actually needed :-( */
  170. /* HBB 970427: seems to have been an incorrect analysis of that error.
  171.  * the problematic plot works without this as well. */
  172. #ifndef OLD_SPLIT_INPLANE
  173. #define OLD_SPLIT_INPLANE 1
  174. #endif
  175.  
  176. /* HBB 970618: this way, new edges introduced by splitting a polygon
  177.  * by the plane of another one will be made visible. Not too useful on
  178.  * production output, but can help in debugging. */
  179. #ifndef SHOW_SPLITTING_EDGES
  180. #define SHOW_SPLITTING_EDGES 0
  181. #endif
  182.  
  183. /********************************/
  184. /* END of configuration section */
  185. /********************************/
  186.  
  187.  
  188.  
  189. /* Variables to hold configurable values. This is meant to prepare for
  190.  * making these settable by the user via 'set hidden [option]...' */
  191.  
  192. int hiddenBacksideLinetypeOffset = BACKSIDE_LINETYPE_OFFSET;
  193. long hiddenTriangleLinesdrawnPattern = TRIANGLE_LINESDRAWN_PATTERN;
  194. int hiddenHandleUndefinedPoints = HANDLE_UNDEFINED_POINTS;
  195. int hiddenShowAlternativeDiagonal = SHOW_ALTERNATIVE_DIAGONAL;
  196. int hiddenHandleBentoverQuadrangles = HANDLE_BENTOVER_QUADRANGLES;
  197.  
  198. /* The functions to map from user 3D space into normalized -1..1 */
  199. #define map_x3d(x) ((x-min_array[FIRST_X_AXIS])*xscale3d-1.0)
  200. #define map_y3d(y) ((y-min_array[FIRST_Y_AXIS])*yscale3d-1.0)
  201. #define map_z3d(z) ((z-floor_z)*zscale3d-1.0)
  202.  
  203. extern int suppressMove;
  204. extern int xright, xleft, ybot, ytop;
  205. extern int xmiddle, ymiddle, xscaler, yscaler;
  206. extern double min_array[], max_array[];
  207. extern int auto_array[], log_array[];
  208. extern double base_array[], log_base_array[];
  209. extern double xscale3d, yscale3d, zscale3d;
  210. extern double floor_z;
  211.  
  212. extern int hidden_no_update, hidden_active;
  213. extern int hidden_line_type_above, hidden_line_type_below;
  214.  
  215. extern double trans_mat[4][4];
  216.  
  217.  
  218. #ifdef HAVE_STRINGIZE
  219. /* ANSI preprocessor concatenation */
  220. # define CONCAT(x,y) x##y
  221. # define CONCAT3(x,y,z) x##y##z
  222. #else
  223. /* K&R preprocessor concatenation */
  224. # define CONCAT(x,y) x/**/y
  225. # define CONCAT3(x,y,z) x/**/y/**/z
  226. #endif
  227.  
  228. /* Bitmap of the screen.  The array for each x value is malloc-ed as needed */
  229. /* HBB 961111: started parametrisation of type t_pnt, to prepare change from
  230.  * short to normal ints in **pnt. The other changes aren't always annotated,
  231.  * so watch out! */
  232. typedef unsigned short int t_pnt;
  233. #define PNT_BITS (CHAR_BIT*sizeof(t_pnt))
  234. /* caution! ANSI-dependant ! ? */
  235. #define PNT_MAX USHRT_MAX
  236. typedef t_pnt *tp_pnt;
  237. static tp_pnt *pnt;
  238.  
  239. /* Testroutine for the bitmap */
  240. /* HBB 961110: fixed slight problem indicated by lclint: 
  241.  * calling IS_UNSET with pnt == 0 (possible?) */
  242. /* HBB 961111: changed for new t_pnt, let compiler take care of optimisation
  243.  * `%' -> `&' */
  244. /* HBB 961124: switched semantics: as this was *always* called as !IS_SET(),
  245.  * it's better to check for *unset* bits right away: */
  246. #define IS_UNSET(X,Y) ((!pnt || pnt[X] == 0) ? 1 : !(((pnt[X])[(Y)/PNT_BITS] >> ((Y)%PNT_BITS)) & 0x01))
  247.  
  248. /* Amount of space we need for one vertical row of bitmap, 
  249.    and byte offset of first used element */
  250. static unsigned long y_malloc;
  251.  
  252. /* These numbers are chosen as dividers into the bitmap. */
  253. static int xfact, yfact;
  254. #define XREDUCE(X) ((X)/xfact)
  255. #define YREDUCE(Y) ((Y)/yfact)
  256.  
  257. /* These variables are only used for drawing the individual polygons
  258.    that make up the 3d figure.  After each polygon is drawn, the
  259.    information is copied to the bitmap: xmin_hl, ymin_hl are used to
  260.    keep track of the range of x values.  The arrays ymin_hl[],
  261.    ymax_hl[] are used to keep track of the minimum and maximum y
  262.    values used for each X value. */
  263.  
  264. /* HBB 961124: new macro, to avoid a wordsize-depence */
  265. #define HL_EXTENT_X_MAX UINT_MAX    /* caution! ANSI-only !? */
  266. static unsigned int xmin_hl, xmax_hl;
  267.  
  268. /* HBB: parametrized type of hl_buffer elements, to allow easier
  269.    changing later on: */
  270. #define HL_EXTENT_Y_MAX SHRT_MAX    /* caution! ANSI-only !? */
  271. typedef short int t_hl_extent_y;
  272. static t_hl_extent_y *ymin_hl, *ymax_hl;
  273.  
  274.  
  275. /* hl_buffer is a buffer which is needed to draw polygons with very small 
  276.    angles correct:
  277.  
  278.    One polygon may be split during the sorting algorithmus into
  279.    several polygons. Before drawing a part of a polygon, I save in
  280.    this buffer all lines of the polygon, which will not yet be drawn.
  281.    If then the actual part draws a virtual line (a invisible line,
  282.    which splits the polygon into 2 parts), it draws the line visible
  283.    in those points, which are set in the buffer.
  284.  
  285.    The layout of the buffer:
  286.    hl_buffer[i] stands for the i'th column of the bitmap.
  287.    Each column is splitted into pairs of points (a,b).
  288.    Each pair describes a cross of one line with the column. */
  289.  
  290. struct Cross {
  291.     int a, b;
  292.     struct Cross GPHUGE *next;
  293. };
  294.  
  295. static struct Cross GPHUGE *GPHUGE * hl_buffer;
  296.  
  297. /* HBB 980303: added a global array of Cross structures, to save lots
  298.  * of gp_alloc() calls (3 millions in a run through all.dem!)  */
  299.  
  300. #define CROSS_STORE_INCREASE 500    /* number of Cross'es to alloc at a time */
  301. static struct Cross *Cross_store = 0;
  302. static int last_Cross_store = 0, free_Cross_store = 0;
  303.  
  304. struct Vertex {
  305.     coordval x, y, z;
  306.     int style;
  307. };
  308.  
  309. /* HBB 971114: two new macros, to properly hide the machinery that
  310.  *implements flagging of vertices as 'undefined' (or 'out of range',
  311.  *handled equally). Thanks to Lars Hecking for pointing this out */
  312. #define FLAG_VERTEX_AS_UNDEFINED(v) \
  313.   do { (v).z = -2.0; } while (0)
  314. #define VERTEX_IS_UNDEFINED(v) ((v).z == -2.0)
  315.  
  316. /* These values serve as flags to detect loops of polygons obscuring
  317.  * each other in in_front() */
  318. typedef enum {
  319.     is_untested = 0, is_tested, is_in_loop, is_moved_or_split
  320. } t_poly_tested;
  321.  
  322. /* Strings for printing out values of type t_poly_tested: */
  323. static const char *string_poly_tested[] =
  324. {
  325.     "is_untested",
  326.     "is_tested",
  327.     "is_in_loop",
  328.     "is_moved_or_split"
  329. };
  330.  
  331. struct Polygon {
  332.     int n;            /* amount of vertices */
  333.     long *vertex;        /* The vertices (indices on vlist) */
  334.     long line;            /* i'th bit shows, if the i'th line should be drawn */
  335.     coordval xmin, xmax, ymin, ymax, zmin, zmax;
  336.     /* min/max in all three directions */
  337.     struct lp_style_type *lp;    /* pointer to l/p properties */
  338.     int style;            /* The style of the lines */
  339.     long next;            /* index of next polygon in z-sorted list */
  340.     long next_frag;        /* next fragment of same polygon... */
  341.     int id;            /* Which polygons belong together? */
  342.     t_poly_tested tested;    /* To determine a loop during the sorting algorithm */
  343.     /* HBB 980317: on the 16 bit PC platforms, the struct size must
  344.      * be a power of two, so it exactly fits into a 64KB segment. First
  345.      * we'll add the TEST_GRIDCHECK fields, regardless wether that
  346.      * feature was activated or not. */
  347. #if TEST_GRIDCHECK || defined(DOS16) || defined(WIN16)
  348.     unsigned int xextent, yextent;
  349. #endif
  350.     /* HBB 980317: the remaining bit of padding. */
  351. #if defined(DOS16) || defined(WIN16) || defined(WIN32)
  352.     char dummies[8];
  353. #endif
  354. };
  355.  
  356. typedef enum {            /* result type for polygon_plane_intersection() */
  357.     infront, inside, behind, intersects
  358. } t_poly_plane_intersect;
  359.  
  360. static struct Vertex GPHUGE *vlist;    /* The vertices */
  361. static struct Polygon GPHUGE *plist;    /* The polygons */
  362. static long nvert, npoly;    /* amount of vertices and polygons */
  363. static long pfree, vert_free;    /* index on the first free entries */
  364.  
  365. static long pfirst;        /* Index on the first polygon */
  366.  
  367. /* macros for (somewhat) safer handling of the dynamic array of
  368.  * polygons, 'plist[]' */
  369. #define EXTEND_PLIST() \
  370.     plist = (struct Polygon GPHUGE *) gp_realloc(plist, \
  371.       (unsigned long)sizeof(struct Polygon)*(npoly += 10L), "hidden plist")
  372.  
  373. #define CHECK_PLIST() if (pfree >= npoly) EXTEND_PLIST()
  374. #define CHECK_PLIST_2() if (pfree+1 >= npoly) EXTEND_PLIST()
  375. #define CHECK_PLIST_EXTRA(extra) if (pfree >= npoly) { EXTEND_PLIST() ; extra; }
  376.  
  377.  
  378. /* precision of calculations in normalized space: */
  379. #define EPSILON 1e-5        /* HBB: was 1.0E-5 */
  380. #define SIGNOF(X)  ( ((X)<-EPSILON) ? -1: ((X)>EPSILON) )
  381.  
  382. /* Some inexact relations: == , > , >= */
  383. #define EQ(X,Y)  (fabs( (X) - (Y) ) < EPSILON)    /* X == Y */
  384. #define GR(X,Y)  ((X) >  (Y) + EPSILON)        /* X >  Y */
  385. #define GE(X,Y)  ((X) >= (Y) - EPSILON)        /* X >= Y */
  386.  
  387. /* Test, if two vertices are near each other */
  388. #define V_EQUAL(a,b) ( GE(0.0, fabs((a)->x - (b)->x) + \
  389.    fabs((a)->y - (b)->y) + fabs((a)->z - (b)->z)) )
  390.  
  391. #define SETBIT(a,i,set) a= (set) ?  (a) | (1L<<(i)) : (a) & ~(1L<<(i))
  392. #define GETBIT(a,i) (((a)>>(i))&1L)
  393.  
  394. #define UINT_BITS (CHAR_BIT*sizeof(unsigned int))
  395. #define USHRT_BITS (CHAR_BIT*sizeof(unsigned short))
  396.  
  397.  
  398. static void print_polygon __PROTO((struct Polygon GPHUGE * p, const char *pname));
  399.  
  400. /* HBB: new routine, to help debug 'Logic errors', mainly */
  401. static void print_polygon(p, pname)
  402. struct Polygon GPHUGE *p;
  403. const char *pname;
  404. {
  405.     struct Vertex GPHUGE *v;
  406.     int i;
  407.  
  408.     fprintf(stderr, "#%s:(ind %d) n:%d, id:%d, next:%ld, tested:%s\n",
  409.         pname, (int) (p - plist), p->n, p->id, p->next,
  410.         string_poly_tested[(int) (p->tested)]);
  411.     fprintf(stderr, "#zmin:%g, zmax:%g\n", p->zmin, p->zmax);
  412.     for (i = 0; i < p->n; i++, v++) {
  413.     v = vlist + p->vertex[i];
  414.     fprintf(stderr, "%g %g %g \t%ld\n", v->x, v->y, v->z, p->vertex[i]);
  415.     }
  416.     /*repeat first vertex, so the line gets closed: */
  417.     v = vlist + p->vertex[0];
  418.     fprintf(stderr, "%g %g %g \t%ld\n", v->x, v->y, v->z, p->vertex[0]);
  419.     /*two blank lines, as a multimesh separator: */
  420.     fputs("\n\n", stderr);
  421. }
  422.  
  423. /* Gets Minimum 'C' value of polygon, C is x, y, or z: */
  424. #define GET_MIN(p, C, min) do { \
  425.   int i; \
  426.   long *v = p->vertex; \
  427.                        \
  428.   min = vlist[*v++].C; \
  429.   for (i = 1; i< p->n; i++, v++) \
  430.     if (vlist[*v].C < min) \
  431.       min = vlist[*v].C; \
  432. } while (0)
  433.  
  434. /* Gets Maximum 'C' value of polygon, C is x, y, or z: */
  435. #define GET_MAX(p, C, max) do { \
  436.   int i; \
  437.   long *v = p->vertex; \
  438.                       \
  439.   max = vlist[*v++].C; \
  440.   for (i = 1; i< p->n; i++, v++) \
  441.     if (vlist[*v].C > max) \
  442.       max = vlist[*v].C; \
  443. } while (0)
  444.  
  445. /* Prototypes for internal functions of this module. */
  446. static void map3d_xyz __PROTO((double x, double y, double z,
  447.                    struct Vertex GPHUGE * v));
  448. static int reduce_polygon __PROTO((int *n, long **v, long *line, int nmin));
  449. static void build_polygon __PROTO((struct Polygon GPHUGE * p, int n,
  450.         long *v, long line, int style, struct lp_style_type * lp,
  451.            long next, long next_frag, int id, t_poly_tested tested));
  452. static GP_INLINE int maybe_build_polygon __PROTO((struct Polygon GPHUGE * p,
  453.      int n, long *v, long line, int style, struct lp_style_type * lp,
  454.            long next, long next_frag, int id, t_poly_tested tested));
  455. static void init_polygons __PROTO((struct surface_points * plots, int pcount));
  456. static int compare_by_zmax __PROTO((const void *p1, const void *p2));
  457. static void sort_by_zmax __PROTO((void));
  458. static int obscured __PROTO((struct Polygon GPHUGE * p));
  459. static GP_INLINE int xy_overlap __PROTO((struct Polygon GPHUGE * a, struct Polygon GPHUGE * b));
  460. static void get_plane __PROTO((struct Polygon GPHUGE * p, double *a, double *b,
  461.                    double *c, double *d));
  462. static t_poly_plane_intersect polygon_plane_intersection __PROTO((struct Polygon GPHUGE * p, double a, double b, double c, double d));
  463. static int intersect __PROTO((struct Vertex GPHUGE * v1,
  464.                   struct Vertex GPHUGE * v2, struct Vertex GPHUGE * v3, struct Vertex GPHUGE * v4));
  465. static int v_intersect __PROTO((struct Vertex GPHUGE * v1,
  466.                 struct Vertex GPHUGE * v2, struct Vertex GPHUGE * v3, struct Vertex GPHUGE * v4));
  467. static int intersect_polygon __PROTO((struct Vertex GPHUGE * v1, struct Vertex GPHUGE * v2, struct Polygon GPHUGE * p));
  468. static int full_xy_overlap __PROTO((struct Polygon GPHUGE * a, struct Polygon GPHUGE * b));
  469. static long build_new_vertex __PROTO((long V1, long V2, double w));
  470. static long split_polygon_by_plane __PROTO((long P, double a, double b,
  471.                     double c, double d, TBOOLEAN f));
  472. static int in_front __PROTO((long Last, long Test));
  473.  
  474. /* HBB 980303: new, for the new back-buffer for *Cross structures: */
  475. static GP_INLINE struct Cross *get_Cross_from_store __PROTO((void));
  476. static GP_INLINE void init_Cross_store __PROTO((void));
  477.  
  478. static GP_INLINE TBOOLEAN hl_buffer_set __PROTO((int xv, int yv));
  479. static GP_INLINE void update_hl_buffer_column __PROTO((int xv, int ya, int yv));
  480. static void draw_clip_line_buffer __PROTO((int x1, int y1, int x2, int y2));
  481. static void draw_clip_line_update __PROTO((int x1, int y1, int x2, int y2,
  482.                        int virtual));
  483. static GP_INLINE void clip_vector_h __PROTO((int x, int y));
  484. static GP_INLINE void clip_vector_virtual __PROTO((int x, int y));
  485. static GP_INLINE void clip_vector_buffer __PROTO((int x, int y));
  486. static void draw __PROTO((struct Polygon GPHUGE * p));
  487.  
  488.  
  489. /* Set the options for hidden3d. To be called from set.c, when the
  490.  * user has begun a command with 'set hidden3d', to parse the rest of
  491.  * that command */
  492. void set_hidden3doptions()
  493. {
  494.     struct value t;
  495.  
  496.     c_token++;
  497.  
  498.     if (END_OF_COMMAND) {
  499.     return;
  500.     }
  501.     if (almost_equals(c_token, "def$aults")) {
  502.     /* reset all parameters to defaults */
  503.     hiddenBacksideLinetypeOffset = BACKSIDE_LINETYPE_OFFSET;
  504.     hiddenTriangleLinesdrawnPattern = TRIANGLE_LINESDRAWN_PATTERN;
  505.     hiddenHandleUndefinedPoints = HANDLE_UNDEFINED_POINTS;
  506.     hiddenShowAlternativeDiagonal = SHOW_ALTERNATIVE_DIAGONAL;
  507.     hiddenHandleBentoverQuadrangles = HANDLE_BENTOVER_QUADRANGLES;
  508.     c_token++;
  509.     if (!END_OF_COMMAND)
  510.         int_error("No further options allowed after 'defaults'", c_token);
  511.     return;
  512.     }
  513.     if (almost_equals(c_token, "off$set")) {
  514.     c_token++;
  515.     hiddenBacksideLinetypeOffset = (int) real(const_express(&t));
  516.     } else if (almost_equals(c_token, "nooff$set")) {
  517.     c_token++;
  518.     hiddenBacksideLinetypeOffset = 0;
  519.     }
  520.     if (END_OF_COMMAND) {
  521.     return;
  522.     }
  523.     if (almost_equals(c_token, "tri$anglepattern")) {
  524.     c_token++;
  525.     hiddenTriangleLinesdrawnPattern = (int) real(const_express(&t));
  526.     }
  527.     if (END_OF_COMMAND) {
  528.     return;
  529.     }
  530.     if (almost_equals(c_token, "undef$ined")) {
  531.     int tmp;
  532.  
  533.     c_token++;
  534.     tmp = (int) real(const_express(&t));
  535.     if (tmp <= 0 || tmp > UNHANDLED)
  536.         tmp = UNHANDLED;
  537.     hiddenHandleUndefinedPoints = tmp;
  538.     } else if (almost_equals(c_token, "nound$efined")) {
  539.     c_token++;
  540.     hiddenHandleUndefinedPoints = UNHANDLED;
  541.     }
  542.     if (END_OF_COMMAND) {
  543.     return;
  544.     }
  545.     if (almost_equals(c_token, "alt$diagonal")) {
  546.     c_token++;
  547.     hiddenShowAlternativeDiagonal = 1;
  548.     } else if (almost_equals(c_token, "noalt$diagonal")) {
  549.     c_token++;
  550.     hiddenShowAlternativeDiagonal = 0;
  551.     }
  552.     if (END_OF_COMMAND) {
  553.     return;
  554.     }
  555.     if (almost_equals(c_token, "bent$over")) {
  556.     c_token++;
  557.     hiddenHandleBentoverQuadrangles = 1;
  558.     } else if (almost_equals(c_token, "nobent$over")) {
  559.     c_token++;
  560.     hiddenHandleBentoverQuadrangles = 0;
  561.     }
  562.     if (!END_OF_COMMAND) {
  563.     int_error("No such option to hidden3d (or wrong order)", c_token);
  564.     }
  565. }
  566.  
  567. /* HBB 970619: new routine for displaying the option settings */
  568. void show_hidden3doptions()
  569. {
  570.     fprintf(stderr,"\
  571. \t  Back side of surfaces has linestyle offset of %d\n\
  572. \t  Bit-Mask of Lines to draw in each triangle is %ld\n\
  573. \t  %d: ",
  574.         hiddenBacksideLinetypeOffset, hiddenTriangleLinesdrawnPattern,
  575.         hiddenHandleUndefinedPoints);
  576.  
  577.     switch (hiddenHandleUndefinedPoints) {
  578.     case OUTRANGE:
  579.     fputs("Outranged and undefined datapoints are omitted from the surface.\n", stderr);
  580.     break;
  581.     case UNDEFINED:
  582.     fputs("Only undefined datapoints are omitted from the surface.\n", stderr);
  583.     break;
  584.     case UNHANDLED:
  585.     fputs("Will not check for undefined datapoints (may cause crashes).\n", stderr);
  586.     break;
  587.     default:
  588.     fputs("This value is illegal!!!\n", stderr);
  589.     break;
  590.     }
  591.     fprintf(stderr,"\
  592. \t  Will %suse other diagonal if it gives a less jaggy outline\n\
  593. \t  Will %sdraw diagonal visibly if quadrangle is 'bent over'\n",
  594.         hiddenShowAlternativeDiagonal ? "" : "not ",
  595.         hiddenHandleBentoverQuadrangles ? "" : "not ");
  596. }
  597.  
  598. /* HBB 971117: new function. */
  599. /* Implements proper 'save'ing of the new hidden3d options... */
  600. void save_hidden3doptions(fp)
  601. FILE *fp;
  602. {
  603.     if (!hidden3d) {
  604.     fputs("set nohidden3d\n", fp);
  605.     return;
  606.     }
  607.     fprintf(fp, "set hidden3d offset %d trianglepattern %ld undefined %d %saltdiagonal %sbentover\n",
  608.         hiddenBacksideLinetypeOffset,
  609.         hiddenTriangleLinesdrawnPattern,
  610.         hiddenHandleUndefinedPoints,
  611.         hiddenShowAlternativeDiagonal ? "" : "no",
  612.         hiddenHandleBentoverQuadrangles ? "" : "no");
  613. }
  614.  
  615. /* Initialize the necessary steps for hidden line removal and
  616.    initialize global variables. */
  617. void init_hidden_line_removal()
  618. {
  619.     int i;
  620.  
  621.     /* Check for some necessary conditions to be set elsewhere: */
  622.     /* HandleUndefinedPoints mechanism depends on these: */
  623.     assert(OUTRANGE == 1);
  624.     assert(UNDEFINED == 2);
  625.     /* '3' makes the test easier in the critical section */
  626.     if (hiddenHandleUndefinedPoints < OUTRANGE)
  627.     hiddenHandleUndefinedPoints = UNHANDLED;
  628.  
  629.     /*  We want to keep the bitmap size less than 2048x2048, so we choose
  630.      *  integer dividers for the x and y coordinates to keep the x and y
  631.      *  ranges less than 2048.  In practice, the x and y sizes for the bitmap
  632.      *  will be somewhere between 1024 and 2048, except in cases where the
  633.      *  coordinates ranges for the device are already less than 1024.
  634.      *  We do this mainly to control the size of the bitmap, but it also
  635.      *  speeds up the computation.  We maintain separate dividers for
  636.      *  x and y.
  637.      */
  638.     xfact = (xright - xleft) / 1024;
  639.     yfact = (ytop - ybot) / 1024;
  640.     if (xfact == 0)
  641.     xfact = 1;
  642.     if (yfact == 0)
  643.     yfact = 1;
  644.     if (pnt == 0) {
  645.     i = XREDUCE(xright) - XREDUCE(xleft) + 1;
  646.     pnt = (tp_pnt *) gp_alloc(
  647.              (unsigned long) (i * sizeof(tp_pnt)), "hidden pnt");
  648.     while (--i >= 0)
  649.         pnt[i] = (tp_pnt) 0;
  650.     }
  651. }
  652.  
  653. /* Reset the hidden line data to a fresh start.                              */
  654. void reset_hidden_line_removal()
  655. {
  656.     int i;
  657.     if (pnt) {
  658.     for (i = 0; i <= XREDUCE(xright) - XREDUCE(xleft); i++) {
  659.         if (pnt[i]) {
  660.         free(pnt[i]);
  661.         pnt[i] = 0;
  662.         };
  663.     };
  664.     };
  665. }
  666.  
  667.  
  668. /* Terminates the hidden line removal process.                  */
  669. /* Free any memory allocated by init_hidden_line_removal above. */
  670. void term_hidden_line_removal()
  671. {
  672.     if (pnt) {
  673.     int j;
  674.     for (j = 0; j <= XREDUCE(xright) - XREDUCE(xleft); j++) {
  675.         if (pnt[j]) {
  676.         free(pnt[j]);
  677.         pnt[j] = 0;
  678.         };
  679.     };
  680.     free(pnt);
  681.     pnt = 0;
  682.     };
  683.     if (ymin_hl)
  684.     free(ymin_hl), ymin_hl = 0;
  685.     if (ymax_hl)
  686.     free(ymax_hl), ymax_hl = 0;
  687. }
  688.  
  689.  
  690. /* Maps from normalized space to terminal coordinates */
  691. #define TERMCOORD(v,x,y) x = ((int)((v)->x * xscaler)) + xmiddle, \
  692.              y = ((int)((v)->y * yscaler)) + ymiddle;
  693.  
  694.  
  695. static void map3d_xyz(x, y, z, v)
  696. double x, y, z;            /* user coordinates */
  697. struct Vertex GPHUGE *v;    /* the point in normalized space */
  698. {
  699.     int i, j;
  700.     double V[4], Res[4],    /* Homogeneous coords. vectors. */
  701.      w = trans_mat[3][3];
  702.     /* Normalize object space to -1..1 */
  703.     V[0] = map_x3d(x);
  704.     V[1] = map_y3d(y);
  705.     V[2] = map_z3d(z);
  706.     V[3] = 1.0;
  707.     Res[0] = 0;            /*HBB 961110: lclint wanted this... Why? */
  708.     for (i = 0; i < 3; i++) {
  709.     Res[i] = trans_mat[3][i];    /* Initiate it with the weight factor. */
  710.     for (j = 0; j < 3; j++)
  711.         Res[i] += V[j] * trans_mat[j][i];
  712.     }
  713.     for (i = 0; i < 3; i++)
  714.     w += V[i] * trans_mat[i][3];
  715.     if (w == 0)
  716.     w = 1.0e-5;
  717.     v->x = Res[0] / w;
  718.     v->y = Res[1] / w;
  719.     v->z = Res[2] / w;
  720. }
  721.  
  722.  
  723. /* remove unnecessary Vertices from a polygon: */
  724. static int reduce_polygon(n, v, line, nmin)
  725. int *n;                /* number of vertices */
  726. long **v;            /* the vertices */
  727. long *line;            /* information which line should be drawn */
  728. int nmin;            /* if the number reduces under nmin, forget polygon */
  729.      /* Return value 1: the reduced polygon has still more or equal 
  730.         vertices than nmin.
  731.         Return value 0: the polygon was trivial; memory is given free */
  732. {
  733.     int less, i;
  734.     long *w = *v;
  735.     for (less = 0, i = 0; i < *n - 1; i++) {
  736.     if (V_EQUAL(vlist + w[i], vlist + w[i + 1])
  737.     /* HBB 970321: try to remove an endless loop bug (part1/2) */
  738.         && ((w[i] == w[i + 1])
  739.         || !GETBIT(*line, i)
  740.         || GETBIT(*line, i + 1)
  741.         || ((i > 0) ? GETBIT(*line, i - 1) : GETBIT(*line, *n - 1))))
  742.         less++;
  743.     else if (less) {
  744.         w[i - less] = w[i];
  745.         SETBIT(*line, i - less, GETBIT(*line, i));
  746.     }
  747.     }
  748.     if (i - less > 0
  749.     && V_EQUAL(vlist + w[i], vlist + w[0])
  750.     /* HBB 970321: try to remove an endless loop bug (part2/2) */
  751.     && (w[i] == w[0]
  752.         || !GETBIT(*line, i)
  753.         || GETBIT(*line, i - 1)
  754.         || GETBIT(*line, 0)))
  755.     less++;
  756.     else if (less) {
  757.     w[i - less] = w[i];
  758.     SETBIT(*line, i - less, GETBIT(*line, i));
  759.     }
  760.     *n -= less;
  761.     if (*n < nmin) {
  762.     free(w);
  763.     *v = 0;            /* HBB 980213: signal that v(==w) is free()d */
  764.     return 0;
  765.     }
  766.     if (less) {
  767.     w = (long *) gp_realloc(w, (unsigned long) sizeof(long) * (*n),
  768.                 "hidden red_poly");
  769.     *v = w;
  770.     for (; less > 0; less--)
  771.         SETBIT(*line, *n + less - 1, 0);
  772.     }
  773.     return 1;
  774. }
  775.  
  776. /* Write polygon in plist */
  777. /*
  778.  * HBB: changed this to precalculate {x|y}{min|max}
  779.  */
  780. static void build_polygon(p, n, v, line, style, lp, next, next_frag, id, tested)
  781. struct Polygon GPHUGE *p;    /* this should point at a free entry in plist */
  782. int n;                /* number of vertices */
  783. long *v;            /* array of indices on the vertices (in vlist) */
  784. long line;            /* information, which line should be drawn */
  785. int style;            /* color of the lines */
  786. struct lp_style_type *lp;    /* pointer to line/point properties */
  787. long next;            /* next polygon in the list */
  788. long next_frag;            /* next fragment of same polygon */
  789. int id;                /* Which polygons belong together? */
  790. t_poly_tested tested;        /* Is polygon already on the right place of list? */
  791. {
  792.     int i;
  793.     struct Vertex vert;
  794.  
  795.     if (p >= plist + npoly)
  796.     fputs("uh oh !\n", stderr);
  797.  
  798.     CHECK_POINTER(plist, p);
  799.  
  800.     p->n = n;
  801.     p->vertex = v;
  802.     p->line = line;
  803.     p->lp = lp;            /* line and point properties */
  804.     p->style = style;
  805.     CHECK_POINTER(vlist, (vlist + v[0]));
  806.     vert = vlist[v[0]];
  807.     p->xmin = p->xmax = vert.x;
  808.     p->ymin = p->ymax = vert.y;
  809.     p->zmin = p->zmax = vert.z;
  810.     for (i = 1; i < n; i++) {
  811.     CHECK_POINTER(vlist, (vlist + v[i]));
  812.     vert = vlist[v[i]];
  813.     if (vert.x < p->xmin)
  814.         p->xmin = vert.x;
  815.     else if (vert.x > p->xmax)
  816.         p->xmax = vert.x;
  817.     if (vert.y < p->ymin)
  818.         p->ymin = vert.y;
  819.     else if (vert.y > p->ymax)
  820.         p->ymax = vert.y;
  821.     if (vert.z < p->zmin)
  822.         p->zmin = vert.z;
  823.     else if (vert.z > p->zmax)
  824.         p->zmax = vert.z;
  825.     }
  826.  
  827. #if TEST_GRIDCHECK
  828. /* FIXME: this should probably use a table of bit-patterns instead... */
  829.     p->xextent = (~(~0U << (unsigned int) ((p->xmax + 1.0) / 2.0 * UINT_BITS + 1.0)))
  830.     & (~0U << (unsigned int) ((p->xmin + 1.0) / 2.0 * UINT_BITS));
  831.     p->yextent = (~(~0U << (unsigned int) ((p->ymax + 1.0) / 2.0 * UINT_BITS + 1.0)))
  832.     & (~0U << (unsigned int) ((p->ymin + 1.0) / 2.0 * UINT_BITS));
  833. #endif
  834.  
  835.     p->next = next;
  836.     p->next_frag = next_frag;    /* HBB 961201: link fragments, to speed up draw() q-loop */
  837.     p->id = id;
  838.     p->tested = tested;
  839. }
  840.  
  841.  
  842. /* get the amount of curves in a plot */
  843. #define GETNCRV(NCRV) do {\
  844.    if(this_plot->plot_type == FUNC3D)        \
  845.      for(icrvs = this_plot->iso_crvs,NCRV = 0; \
  846.      icrvs;icrvs = icrvs->next,NCRV++);  \
  847.    else if(this_plot->plot_type == DATA3D) \
  848.         NCRV = this_plot->num_iso_read;    \
  849.     else {                                 \
  850.         graph_error("Plot type is neither function nor data"); \
  851.         return; /* stops a gcc -Wunitialised warning */        \
  852.     } \
  853. } while (0)
  854.  
  855. /* Do we see the top or bottom of the polygon, or is it 'on edge'? */
  856. #define GET_SIDE(vlst,csign) do { \
  857.   double c = vlist[vlst[0]].x * (vlist[vlst[1]].y - vlist[vlst[2]].y) + \
  858.     vlist[vlst[1]].x * (vlist[vlst[2]].y - vlist[vlst[0]].y) + \
  859.     vlist[vlst[2]].x * (vlist[vlst[0]].y - vlist[vlst[1]].y); \
  860.   csign = SIGNOF (c); \
  861. } while (0)
  862.  
  863. /* HBB 970131: new function, to allow handling of undefined datapoints
  864.  * by leaving out the corresponding polygons from the list.
  865.  * Return Value: 1 if polygon was built, 0 otherwise */
  866. static GP_INLINE int maybe_build_polygon(p, n, v, line, style, lp,
  867.                      next, next_frag, id, tested)
  868. struct Polygon GPHUGE *p;
  869. struct lp_style_type *lp;
  870. int n, style, id;
  871. t_poly_tested tested;
  872. long *v, line, next, next_frag;
  873. {
  874.     int i;
  875.  
  876.     assert(v);
  877.  
  878.     if (hiddenHandleUndefinedPoints != UNHANDLED)
  879.     for (i = 0; i < n; i++)
  880.         if (VERTEX_IS_UNDEFINED(vlist[v[i]])) {
  881.         free(v);    /* HBB 980213: free mem... */
  882.         v = 0;
  883.         return 0;    /* *don't* build the polygon! */
  884.         }
  885.     build_polygon(p, n, v, line, style, lp, next, next_frag, id, tested);
  886.     return 1;
  887. }
  888.  
  889. /* HBB 970322: new macro, invented because code like this occured several
  890.  * times inside init_polygons, esp. after I added the option to choose the
  891.  * other diagonal to divide a quad when necessary */
  892. /* HBB 980213: Here and below in init_polygons, added code to properly
  893.  * free any allocated vertex lists ('v' here), or avoid allocating
  894.  * them in the first place. Thanks to Stefan Schroepfer for reporting
  895.  * the leak.  */
  896. #define PREPARE_POLYGON(n,v,i0,i1,i2,line,c,border,i_chk,ncrv_chk,style) do {\
  897.   (n) = 3; \
  898.   if (VERTEX_IS_UNDEFINED(vlist[vert_free + (i0)])\
  899.       || VERTEX_IS_UNDEFINE-”-,-X-e + (i1)]) \
  900.       || VERTEX_IS_UNDEFINED(vlist[vert_free + (i2)])) { \
  901.     /* These values avoid any further action on this polygon */\
  902.     (v) = 0; /* flag this as undefined */ \
  903.     (c) = (border) = 0; \
  904.   } else {\
  905.   (v) = (long *) gp_alloc ((unsigned long) sizeof (long) * (n), "hidden PREPARE_POLYGON"); \
  906.   (v)[0] = vert_free + (i0);\
  907.   (v)[1] = vert_free + (i1);\
  908.   (v)[2] = vert_free + (i2);\
  909.   (line) = hiddenTriangleLinesdrawnPattern;\
  910.     GET_SIDE((v),(c));\
  911.     /* Is this polygon at the border of the surface? */\
  912.     (border) = (i == (i_chk) || ncrv == (ncrv_chk));\
  913.   }\
  914.   (style) = ((c) >= 0) ? above : below;\
  915. } while (0)
  916.  
  917. static void init_polygons(plots, pcount)
  918. struct surface_points *plots;
  919. int pcount;
  920.      /* Initialize vlist, plist */
  921. {
  922.     long int i;
  923.     struct surface_points *this_plot;
  924.     int surface;
  925.     long int ncrv, ncrv1;
  926.     struct iso_curve *icrvs;
  927.     int row_offset;
  928.     int id;
  929.     int n1, n2;            /* number of vertices of a Polygon */
  930.     long *v1, *v2;        /* the Vertices */
  931.     long line1, line2;        /* which  lines should be drawn */
  932.     int above = -1, below, style1, style2;    /* the line color */
  933.     struct lp_style_type *lp;    /* pointer to line and point properties */
  934.     int c1, c2;            /* do we see the front or the back */
  935.     TBOOLEAN border1, border2;    /* is it at the border of the surface */
  936.  
  937.     /* allocate memory for polylist and nodelist */
  938.     nvert = npoly = 0;
  939.     for (this_plot = plots, surface = 0;
  940.      surface < pcount;
  941.      this_plot = this_plot->next_sp, surface++) {
  942.     GETNCRV(ncrv);
  943.     switch (this_plot->plot_style) {
  944.     case LINESPOINTS:
  945.     case STEPS:        /* handle as LINES */
  946.     case FSTEPS:
  947.     case HISTEPS:
  948.     case LINES:
  949. #if FOUR_TRIANGLES
  950.         nvert += ncrv * (2 * this_plot->iso_crvs->p_count - 1);
  951. #else
  952.         nvert += ncrv * (this_plot->iso_crvs->p_count);
  953. #endif
  954.         npoly += 2 * (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  955.         break;
  956.     case DOTS:
  957.     case XERRORBARS:    /* handle as POINTSTYLE */
  958.     case YERRORBARS:
  959.     case XYERRORBARS:
  960.     case BOXXYERROR:
  961.     case BOXERROR:
  962.     case CANDLESTICKS:    /* HBB: these as well */
  963.     case FINANCEBARS:
  964.     case VECTOR:
  965.     case POINTSTYLE:
  966.         nvert += ncrv * (this_plot->iso_crvs->p_count);
  967.         npoly += (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  968.         break;
  969.     case BOXES:        /* handle as IMPULSES */
  970.     case IMPULSES:
  971.         nvert += 2 * ncrv * (this_plot->iso_crvs->p_count);
  972.         npoly += (ncrv - 1) * (this_plot->iso_crvs->p_count - 1);
  973.         break;
  974.     }
  975.     }
  976.     vlist = (struct Vertex GPHUGE *)
  977.     gp_alloc((unsigned long) sizeof(struct Vertex) * nvert, "hidden vlist");
  978.     plist = (struct Polygon GPHUGE *)
  979.     gp_alloc((unsigned long) sizeof(struct Polygon) * npoly, "hidden plist");
  980.  
  981.     /* initialize vlist: */
  982.     for (vert_free = 0, this_plot = plots, surface = 0;
  983.      surface < pcount;
  984.      this_plot = this_plot->next_sp, surface++) {
  985.     switch (this_plot->plot_style) {
  986.     case LINESPOINTS:
  987.     case BOXERROR:        /* handle as POINTSTYLE */
  988.     case BOXXYERROR:
  989.     case XERRORBARS:
  990.     case YERRORBARS:
  991.     case XYERRORBARS:
  992.     case CANDLESTICKS:    /* HBB: these as well */
  993.     case FINANCEBARS:
  994.     case VECTOR:
  995.     case POINTSTYLE:
  996.         above = this_plot->lp_properties.p_type;
  997.         break;
  998.     case STEPS:        /* handle as LINES */
  999.     case FSTEPS:
  1000.     case HISTEPS:
  1001.     case LINES:
  1002.     case DOTS:
  1003.     case BOXES:        /* handle as IMPULSES */
  1004.     case IMPULSES:
  1005.         above = -1;
  1006.         break;
  1007.     }
  1008.     GETNCRV(ncrv1);
  1009.     for (ncrv = 0, icrvs = this_plot->iso_crvs;
  1010.          ncrv < ncrv1 && icrvs;
  1011.          ncrv++, icrvs = icrvs->next) {
  1012.         struct coordinate GPHUGE *points = icrvs->points;
  1013.  
  1014.         for (i = 0; i < icrvs->p_count; i++) {
  1015.         /* As long as the point types keep OUTRANGE == 1 and
  1016.          * UNDEFINED == 2, we can just compare
  1017.          * hiddenHandleUndefinedPoints to the type of the point. To
  1018.          * simplify this, let UNHANDLED := UNDEFINED+1 mean 'do not
  1019.          * handle undefined or outranged points at all' (dangerous
  1020.          * option, that is...) */
  1021.         if (points[i].type >= hiddenHandleUndefinedPoints) {
  1022.             /* mark this vertex as a bad one */
  1023.             FLAG_VERTEX_AS_UNDEFINED(vlist[vert_free++]);
  1024.             continue;
  1025.         }
  1026.         map3d_xyz(points[i].x, points[i].y, points[i].z, vlist + vert_free);
  1027.         vlist[vert_free++].style = above;
  1028.         if (this_plot->plot_style == IMPULSES
  1029.             || this_plot->plot_style == BOXES) {
  1030.             map3d_xyz(points[i].x, points[i].y, min_array[FIRST_Z_AXIS], vlist + vert_free);
  1031.             vlist[vert_free++].style = above;
  1032.         }
  1033. #if FOUR_TRIANGLES
  1034.         /* FIXME: handle other styles as well! */
  1035.         if (this_plot->plot_style == LINES &&
  1036.             i < icrvs->p_count - 1)
  1037.             vert_free++;    /* keep one entry free for quad-center */
  1038. #endif
  1039.         }
  1040.     }
  1041.     }
  1042.  
  1043.     /* initialize plist: */
  1044.     id = 0;
  1045.     for (pfree = vert_free = 0, this_plot = plots, surface = 0;
  1046.      surface < pcount;
  1047.      this_plot = this_plot->next_sp, surface++) {
  1048.     row_offset = this_plot->iso_crvs->p_count;
  1049.     lp = &(this_plot->lp_properties);
  1050.     above = this_plot->lp_properties.l_type;
  1051.     below = this_plot->lp_properties.l_type + hiddenBacksideLinetypeOffset;
  1052.     GETNCRV(ncrv1);
  1053.     for (ncrv = 0, icrvs = this_plot->iso_crvs;
  1054.          ncrv < ncrv1 && icrvs;
  1055.          ncrv++, icrvs = icrvs->next)
  1056.         for (i = 0; i < icrvs->p_count; i++)
  1057.         switch (this_plot->plot_style) {
  1058.         case LINESPOINTS:
  1059.         case STEPS:    /* handle as LINES */
  1060.         case FSTEPS:
  1061.         case HISTEPS:
  1062.         case LINES:
  1063.             if (i < icrvs->p_count - 1 && ncrv < ncrv1 - 1) {
  1064. #if FOUR_TRIANGLES
  1065.             long *v3, *v4;
  1066.             int n3, n4;
  1067.             long line3, line4;
  1068.             int c3, c4, style3, style4;
  1069.             TBOOLEAN border3, border4, inc_id = FALSE;
  1070.  
  1071.             /* Build medium vertex: */
  1072.             vlist[vert_free + 1].x =
  1073.                 (vlist[vert_free].x + vlist[vert_free + 2 * row_offset - 1].x +
  1074.                  vlist[vert_free + 2].x + vlist[vert_free + 2 * row_offset + 1].x
  1075.                 ) / 4;
  1076.             vlist[vert_free + 1].y =
  1077.                 (vlist[vert_free].y + vlist[vert_free + 2 * row_offset - 1].y +
  1078.                  vlist[vert_free + 2].y + vlist[vert_free + 2 * row_offset + 1].y
  1079.                 ) / 4;
  1080.             vlist[vert_free + 1].z =
  1081.                 (vlist[vert_free].z + vlist[vert_free + 2 * row_offset - 1].z +
  1082.                  vlist[vert_free + 2].z + vlist[vert_free + 2 * row_offset + 1].z
  1083.                 ) / 4;
  1084.             vlist[vert_free + 1].style = above;
  1085.  
  1086.             PREPARE_POLYGON(n1, v1, 2, 0, 1, line1, c1,
  1087.                     border1, 0, -1, style1);
  1088.             /* !!FIXME!! special casing (see the older code) still needs to be
  1089.              * implemented! */
  1090.             if ((c1 || border1)
  1091.                 && reduce_polygon(&n1, &v1, &line1, 3)) {
  1092.                 CHECK_PLIST();
  1093.                 pfree += maybe_build_polygon(plist + pfree, n1, v1, line1,
  1094.                           style1, lp, -1L, -1L, id++,
  1095.                              is_untested);
  1096.                 inc_id = TRUE;
  1097.             }
  1098.             PREPARE_POLYGON(n2, v2, 0, 2 * row_offset - 1, 1, line2, c2,
  1099.                     border2, -1, 0, style2);
  1100.             if ((c2 || border2)
  1101.                 && reduce_polygon(&n2, &v2, &line2, 3)) {
  1102.                 CHECK_PLIST();
  1103.                 pfree += maybe_build_polygon(plist + pfree, n2, v2, line2,
  1104.                           style2, lp, -1L, -1L, id++,
  1105.                              is_untested);
  1106.                 inc_id = TRUE;
  1107.             }
  1108.             PREPARE_POLYGON(n3, v3, 2 * row_offset - 1, 2 * row_offset + 1, 1,
  1109.                   line3, c3, border3, icrvs->p_count - 2, -1,
  1110.                     style3);
  1111.             if ((c3 || border3)
  1112.                 && reduce_polygon(&n3, &v3, &line3, 3)) {
  1113.                 CHECK_PLIST();
  1114.                 pfree += maybe_build_polygon(plist + pfree, n3, v3, line3,
  1115.                           style3, lp, -1L, -1L, id++,
  1116.                              is_untested);
  1117.                 inc_id = TRUE;
  1118.             }
  1119.             PREPARE_POLYGON(n4, v4, 2 * row_offset + 1, 0, 1, line4, c4,
  1120.                     border4, -1, ncrv1 - 2, style4);
  1121.             if ((c4 || border4)
  1122.                 && reduce_polygon(&n4, &v4, &line4, 3)) {
  1123.                 CHECK_PLIST();
  1124.                 pfree += maybe_build_polygon(plist + pfree, n4, v4, line4,
  1125.                           style4, lp, -1L, -1L, id++,
  1126.                              is_untested);
  1127.                 inc_id = TRUE;
  1128.             }
  1129.             /*if (inc_id)
  1130.                id++; */
  1131.             vert_free++;
  1132.  
  1133. #else /* FOUR_TRIANGLES */
  1134.  
  1135.             PREPARE_POLYGON(n1, v1, row_offset, 0, 1, line1, c1,
  1136.                     border1, 0, 0, style1);
  1137.             PREPARE_POLYGON(n2, v2, 1, row_offset + 1, row_offset, line2,
  1138.                   c2, border2, icrvs->p_count - 2, ncrv1 - 2,
  1139.                     style2);
  1140.  
  1141.             /* if this makes at least one of the two triangles
  1142.              * visible, use the other diagonal here */
  1143.             if (hiddenShowAlternativeDiagonal
  1144.                 && !(v1)    /* HBB 980213: only try this in the case of */
  1145.                 &&!(v2)    /* undefined vertices -> *missing* trigs */
  1146.                 ) {
  1147.                 if (VERTEX_IS_UNDEFINED(vlist[vert_free + 1])) {
  1148.                 PREPARE_POLYGON(n1, v1, row_offset + 1, row_offset, 0,
  1149.                     line1, c1, border1, 0, ncrv1 - 2,
  1150.                         style1);
  1151.                 } else if (VERTEX_IS_UNDEFINED(vlist[vert_free + row_offset])) {
  1152.                 PREPARE_POLYGON(n1, v1, 0, 1, row_offset + 1, line1,
  1153.                       c1, border1, icrvs->p_count - 2, 0,
  1154.                         style1);
  1155.                 }
  1156.             }
  1157.             /* HBB 970323: if visible sides of the two triangles
  1158.              * differs, make diagonal visible! */
  1159.             /* HBB 970428: FIXME: should this be changed to only
  1160.              * activate *one* of the triangles' diagonals? If so, how
  1161.              * to find out the right one? Maybe the one with its
  1162.              * normal closer to the z direction? */
  1163.             if (hiddenHandleBentoverQuadrangles
  1164.                 && (c1 * c2) < 0
  1165.                 ) {
  1166.                 SETBIT(line1, n1 - 1, 1);
  1167.                 SETBIT(line2, n2 - 1, 1);
  1168.             }
  1169.             if ((c1 || border1)
  1170.                 && reduce_polygon(&n1, &v1, &line1, 3)) {
  1171.                 if ((c2 || border2)
  1172.                 && reduce_polygon(&n2, &v2, &line2, 3)) {
  1173.                 /* These two will need special handling, to
  1174.                  * ensure the links from one to the other
  1175.                  * are only set up when *both* polygons are
  1176.                  * valid */
  1177.                 int r1, r2;    /* temp. results */
  1178.  
  1179.                 CHECK_PLIST_2();
  1180.                 r1 = maybe_build_polygon(plist + pfree, n1, v1,
  1181.                           line1, style1, lp, -1L,
  1182.                            -1L, id, is_untested);
  1183.                 r2 = maybe_build_polygon(plist + pfree + r1, n2, v2,
  1184.                           line2, style2, lp, -1L,
  1185.                          -1L, id++, is_untested);
  1186.                 if (r1 && r2) {
  1187.                     plist[pfree].next_frag = pfree + 1;
  1188.                     plist[pfree + 1].next_frag = pfree;
  1189.                 }
  1190.                 pfree += r1 + r2;
  1191.                 } else {    /* P1 was ok, but P2 wasn't */
  1192. /* HBB 980213: P2 is invisible, so free its vertex list: */
  1193.                 if (v2) {
  1194.                     free(v2);
  1195.                     v2 = 0;
  1196.                 }
  1197. /* HBB 970323: if other polygon wasn't drawn, draw this polygon's
  1198.  * diagonal visibly (not only if it was 'on edge'). */
  1199.                 /*if (c2 || border2) */
  1200.                 SETBIT(line1, n1 - 1, 1);
  1201.                 CHECK_PLIST();
  1202.                 pfree += maybe_build_polygon(plist + pfree, n1, v1, line1,
  1203.                           style1, lp, -1L, -1L, id++,
  1204.                                  is_untested);
  1205.                 }
  1206.             } else {    /* P1 was not ok */
  1207. /* HBB 980213: P1 is invisible, so free its vertex list: */
  1208.                 if (v1) {
  1209.                 free(v1);
  1210.                 v1 = 0;
  1211.                 }
  1212. /* HBB 970323: if other polygon wasn't drawn, draw this polygon's
  1213.  * diagonal visibly (not only if it was 'on edge'). */
  1214.                 /*if (c1 || border1) */
  1215.                 SETBIT(line2, n2 - 1, 1);
  1216. /* HBB 970522: yet another problem triggered by the handling of
  1217.  * undefined points: if !(c2||border2), we mustn't try to reduce
  1218.  * triangle 2. */
  1219.                 if (0
  1220.                 || (1
  1221.                     && (c2 || border2)
  1222.                   && reduce_polygon(&n2, &v2, &line2, 2))
  1223.                 || (1
  1224. /* HBB 980213: reduce_... above might have zeroed v... */
  1225.                     && (v2)
  1226.                     && (c1 || border1))) {
  1227.                 CHECK_PLIST();
  1228.                 pfree += maybe_build_polygon(plist + pfree, n2, v2, line2,
  1229.                           style2, lp, -1L, -1L, id++,
  1230.                                  is_untested);
  1231.                 } else {
  1232.                 /* HBB 980213: sigh... both P1 and P2 were invisible, but
  1233.                  * at least one of them was a not undefined.. */
  1234.                 if (v1)
  1235.                     free(v1);
  1236.                 if (v2)
  1237.                     free(v2);
  1238.                 }
  1239.             }
  1240. #endif /* else: FOUR_TRIANGLES */
  1241.             }
  1242.             vert_free++;
  1243.             break;
  1244.         case DOTS:
  1245.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * 1,
  1246.                        "hidden v1 for dots");
  1247.             v1[0] = vert_free++;
  1248.             CHECK_PLIST();
  1249.             pfree += maybe_build_polygon(plist + pfree, 1, v1, 1L, above,
  1250.                     lp, -1L, -1L, id++, is_untested);
  1251.             break;
  1252.         case BOXERROR:    /* handle as POINTSTYLE */
  1253.         case BOXXYERROR:
  1254.         case XERRORBARS:    /* handle as POINTSTYLE */
  1255.         case YERRORBARS:
  1256.         case XYERRORBARS:
  1257.         case CANDLESTICKS:    /* HBB: these as well */
  1258.         case FINANCEBARS:
  1259.         case POINTSTYLE:
  1260.         case VECTOR:
  1261.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * 1,
  1262.                        "hidden v1 for point");
  1263.             v1[0] = vert_free++;
  1264.             CHECK_PLIST();
  1265.             pfree += maybe_build_polygon(plist + pfree, 1, v1, 0L, above,
  1266.                     lp, -1L, -1L, id++, is_untested);
  1267.             break;
  1268.         case BOXES:    /* handle as IMPULSES */
  1269.         case IMPULSES:
  1270.             n1 = 2;
  1271.             v1 = (long *) gp_alloc((unsigned long) sizeof(long) * n1,
  1272.                        "hidden v1 for impulse");
  1273.             v1[0] = vert_free++;
  1274.             v1[1] = vert_free++;
  1275.             line1 = 2L;
  1276.             (void) reduce_polygon(&n1, &v1, &line1, 1);
  1277.             CHECK_PLIST();
  1278.             pfree += maybe_build_polygon(plist + pfree, n1, v1, line1, above,
  1279.                     lp, -1L, -1L, id++, is_untested);
  1280.             break;
  1281.         }
  1282.     }
  1283. }
  1284.  
  1285. static int compare_by_zmax(p1, p2)
  1286. const void *p1, *p2;
  1287. {
  1288.     return (SIGNOF(plist[*(const long *) p2].zmax - plist[*(const long *) p1].zmax));
  1289. }
  1290.  
  1291. static void sort_by_zmax()
  1292.      /* and build the list (return_value = start of list) */
  1293. {
  1294.     long *sortarray, i;
  1295.     struct Polygon GPHUGE *this;
  1296.     sortarray = (long *) gp_alloc((unsigned long) sizeof(long) * pfree, "hidden sortarray");
  1297.     for (i = 0; i < pfree; i++)
  1298.     sortarray[i] = i;
  1299.     qsort(sortarray, (size_t) pfree, sizeof(long), compare_by_zmax);
  1300.     this = plist + sortarray[0];
  1301.     for (i = 1; i < pfree; i++) {
  1302.     this->next = sortarray[i];
  1303.     this = plist + sortarray[i];
  1304.     }
  1305.     this->next = -1L;
  1306.     pfirst = sortarray[0];
  1307.     free(sortarray);
  1308. }
  1309.  
  1310.  
  1311. /* HBB: try if converting this code to unsigned int (from signed shorts)
  1312.  *  fixes any of the remaining bugs. In the same motion, remove
  1313.  *  hardcoded sizeof(short) = 2 (09.10.1996) */
  1314.  
  1315. #define LASTBITS (USHRT_BITS -1)    /* ????? */
  1316. static int obscured(p)
  1317.      /* is p obscured by already drawn polygons? (only a simple minimax-test) */
  1318. struct Polygon GPHUGE *p;
  1319. {
  1320.     int l_xmin, l_xmax, l_ymin, l_ymax;        /* HBB 961110: avoid shadowing external names */
  1321.     t_pnt mask1, mask2;
  1322.     long indx1, indx2, k, m;
  1323.     tp_pnt cpnt;
  1324.     /* build the minimax-box */
  1325.     l_xmin = (p->xmin * xscaler) + xmiddle;
  1326.     l_xmax = (p->xmax * xscaler) + xmiddle;
  1327.     l_ymin = (p->ymin * yscaler) + ymiddle;
  1328.     l_ymax = (p->ymax * yscaler) + ymiddle;
  1329.     if (l_xmin < xleft)
  1330.     l_xmin = xleft;
  1331.     if (l_xmax > xright)
  1332.     l_xmax = xright;
  1333.     if (l_ymin < ybot)
  1334.     l_ymin = ybot;
  1335.     if (l_ymax > ytop)
  1336.     l_ymax = ytop;
  1337.     if (l_xmin > l_xmax || l_ymin > l_ymax)
  1338.     return 1;        /* not viewable */
  1339.     l_ymin = YREDUCE(l_ymin) - YREDUCE(ybot);
  1340.     l_ymax = YREDUCE(l_ymax) - YREDUCE(ybot);
  1341.     l_xmin = XREDUCE(l_xmin) - XREDUCE(xleft);
  1342.     l_xmax = XREDUCE(l_xmax) - XREDUCE(xleft);
  1343.     /* Now check bitmap */
  1344.     indx1 = l_ymin / PNT_BITS;
  1345.     indx2 = l_ymax / PNT_BITS;
  1346.     mask1 = PNT_MAX << (((unsigned) l_ymin) % PNT_BITS);
  1347.     mask2 = PNT_MAX >> ((~((unsigned) l_ymax)) % PNT_BITS);
  1348.     /* HBB: lclint wanted this: */
  1349.     assert(pnt != 0);
  1350.     for (m = l_xmin; m <= l_xmax; m++) {
  1351.     if (pnt[m] == 0)
  1352.         return 0;        /* not obscured */
  1353.     cpnt = pnt[m] + indx1;
  1354.     if (indx1 == indx2) {
  1355.         if (~(*cpnt) & mask1 & mask2)
  1356.         return 0;
  1357.     } else {
  1358.         if (~(*cpnt++) & mask1)
  1359.         return 0;
  1360.         for (k = indx1 + 1; k < indx2; k++)
  1361.         if ((*cpnt++) != PNT_MAX)
  1362.             return 0;
  1363.         if (~(*cpnt++) & mask2)
  1364.         return 0;
  1365.     }
  1366.     }
  1367.     return 1;
  1368. }
  1369.  
  1370.  
  1371. void draw_line_hidden(x1, y1, x2, y2)
  1372. unsigned int x1, y1, x2, y2;
  1373. {
  1374.     register struct termentry *t = term;
  1375.     TBOOLEAN flag;
  1376.     register int xv, yv, errx, erry, err;
  1377.     register unsigned int xvr, yvr;
  1378.     int unsigned xve, yve;
  1379.     register int dy, nstep = 0, dyr;
  1380.     int i;
  1381.  
  1382.     if (x1 > x2) {
  1383.     xvr = x2;
  1384.     yvr = y2;
  1385.     xve = x1;
  1386.     yve = y1;
  1387.     } else {
  1388.     xvr = x1;
  1389.     yvr = y1;
  1390.     xve = x2;
  1391.     yve = y2;
  1392.     };
  1393.     errx = XREDUCE(xve) - XREDUCE(xvr);
  1394.     erry = YREDUCE(yve) - YREDUCE(yvr);
  1395.     dy = (erry > 0 ? 1 : -1);
  1396.     dyr = dy * yfact;
  1397.     switch (dy) {
  1398.     case 1:
  1399.     nstep = errx + erry;
  1400.     errx = -errx;
  1401.     break;
  1402.     case -1:
  1403.     nstep = errx - erry;
  1404.     errx = -errx;
  1405.     erry = -erry;
  1406.     break;
  1407.     };
  1408.     err = errx + erry;
  1409.     errx <<= 1;
  1410.     erry <<= 1;
  1411.     xv = XREDUCE(xvr) - XREDUCE(xleft);
  1412.     yv = YREDUCE(yvr) - YREDUCE(ybot);
  1413.     (*t->move) (xvr, yvr);
  1414.     flag = !IS_UNSET(xv, yv);
  1415.     if (!hidden_no_update) {    /* Check first point */
  1416.     assert(ymax_hl != 0);
  1417.     assert(ymin_hl != 0);
  1418.     if (xv < xmin_hl)
  1419.         xmin_hl = xv;
  1420.     if (xv > xmax_hl)
  1421.         xmax_hl = xv;
  1422.     if (yv > ymax_hl[xv])
  1423.         ymax_hl[xv] = yv;
  1424.     if (yv < ymin_hl[xv])
  1425.         ymin_hl[xv] = yv;
  1426.     };
  1427.     for (i = 0; i < nstep; i++) {
  1428.     if (err < 0) {
  1429.         xv++;
  1430.         xvr += xfact;
  1431.         err += erry;
  1432.     } else {
  1433.         yv += dy;
  1434.         yvr += dyr;
  1435.         err += errx;
  1436.     };
  1437.     if (IS_UNSET(xv, yv)) {
  1438.         if (flag) {
  1439.         (*t->move) (xvr, yvr);
  1440.         flag = FALSE;
  1441.         };
  1442.     } else {
  1443.         if (!flag) {
  1444.         (*t->vector) (xvr, yvr);
  1445.         flag = TRUE;
  1446.         };
  1447.     };
  1448.     if (!hidden_no_update) {
  1449.         /* HBB 961110: lclint wanted these: */
  1450.         assert(ymax_hl != 0);
  1451.         assert(ymin_hl != 0);
  1452.         if (xv < xmin_hl)
  1453.         xmin_hl = xv;
  1454.         if (xv > xmax_hl)
  1455.         xmax_hl = xv;
  1456.         if (yv > ymax_hl[xv])
  1457.         ymax_hl[xv] = yv;
  1458.         if (yv < ymin_hl[xv])
  1459.         ymin_hl[xv] = yv;
  1460.     };
  1461.     };
  1462.     if (!flag)
  1463.     (*t->vector) (xve, yve);
  1464.     return;
  1465. }
  1466.  
  1467.  
  1468. static GP_INLINE int xy_overlap(a, b)
  1469.      /* Do a and b overlap in x or y (minimax test) */
  1470. struct Polygon GPHUGE *a, GPHUGE * b;
  1471. {
  1472. #if TEST_GRIDCHECK
  1473.     /* First, check by comparing the extent bit patterns: */
  1474.     if (!((a->xextent & b->xextent) && (a->yextent & b->yextent)))
  1475.     return 0;
  1476. #endif
  1477.     return ((int) (GE(b->xmax, a->xmin) && GE(a->xmax, b->xmin)
  1478.            && GE(b->ymax, a->ymin) && GE(a->ymax, b->ymin)));
  1479. }
  1480.  
  1481.  
  1482. static void get_plane(p, a, b, c, d)
  1483. struct Polygon GPHUGE *p;
  1484. double *a, *b, *c, *d;
  1485. {
  1486.     int i;
  1487.     struct Vertex GPHUGE *v1, GPHUGE * v2;
  1488.     double x, y, z, s;
  1489.     if (p->n == 1) {
  1490.     *a = 0.0;
  1491.     *b = 0.0;
  1492.     *c = 1.0;
  1493.     *d = -vlist[p->vertex[0]].z;
  1494.     return;
  1495.     }
  1496.     v1 = vlist + p->vertex[p->n - 1];
  1497.     v2 = vlist + p->vertex[0];
  1498.     *a = (v1->y - v2->y) * (v1->z + v2->z);
  1499.     *b = (v1->z - v2->z) * (v1->x + v2->x);
  1500.     *c = (v1->x - v2->x) * (v1->y + v2->y);
  1501.     for (i = 0; i < p->n - 1; i++) {
  1502.     v1 = vlist + p->vertex[i];
  1503.     v2 = vlist + p->vertex[i + 1];
  1504.     *a += (v1->y - v2->y) * (v1->z + v2->z);
  1505.     *b += (v1->z - v2->z) * (v1->x + v2->x);
  1506.     *c += (v1->x - v2->x) * (v1->y + v2->y);
  1507.     }
  1508.     s = sqrt(*a * *a + *b * *b + *c * *c);
  1509.     if (GE(0.0, s)) {        /* => s==0 => the vertices are in one line */
  1510.     v1 = vlist + p->vertex[0];
  1511.     for (i = 1; i < p->n; i++) {
  1512.         v2 = vlist + p->vertex[i];
  1513.         if (!V_EQUAL(v1, v2))
  1514.         break;
  1515.     }
  1516.     /* (x,y,z) is l.u. from <v1, v2> */
  1517.     x = v1->x;
  1518.     y = v1->y;
  1519.     z = v1->z;
  1520.     if (EQ(y, v2->y))
  1521.         y += 1.0;
  1522.     else
  1523.         x += 1.0;
  1524.     /* build a vector that is orthogonal to the line of the polygon */
  1525.     *a = v1->y * (v2->z - z) + v2->y * (z - v1->z) + y * (v1->z - v2->z);
  1526.     *b = v1->z * (v2->x - x) + v2->z * (x - v1->x) + z * (v1->x - v2->x);
  1527.     *c = v1->x * (v2->y - y) + v2->x * (y - v1->y) + x * (v1->y - v2->y);
  1528.     s = sqrt(*a * *a + *b * *b + *c * *c);
  1529.     }
  1530.     if (*c < 0.0)        /* ensure that normalized c is > 0 */
  1531.     s *= -1.0;        /* The exact relation, because c can be very small */
  1532.     *a /= s;
  1533.     *b /= s;
  1534.     *c /= s;
  1535.     *d = -*a * v1->x - *b * v1->y - *c * v1->z;
  1536.     return;
  1537. }
  1538.  
  1539.  
  1540. static t_poly_plane_intersect
  1541.  polygon_plane_intersection(p, a, b, c, d)
  1542. struct Polygon GPHUGE *p;
  1543. double a, b, c, d;
  1544. {
  1545.     int i, sign, max_sign, min_sign;
  1546.     struct Vertex GPHUGE *v;
  1547.  
  1548.     CHECK_POINTER(plist, p);
  1549.  
  1550.     v = vlist + p->vertex[0];
  1551.     max_sign = min_sign = SIGNOF(a * v->x + b * v->y + c * v->z + d);
  1552.     for (i = 1; i < p->n; i++) {
  1553.     v = vlist + p->vertex[i];
  1554.     sign = SIGNOF(a * v->x + b * v->y + c * v->z + d);
  1555.     if (sign > max_sign)
  1556.         max_sign = sign;
  1557.     else if (sign < min_sign)
  1558.         min_sign = sign;
  1559.     }
  1560.  
  1561.     /* Yes, this is a bit tricky, but it's simpler for the computer... */
  1562.     if (min_sign == -1) {
  1563.     if (max_sign == 1)
  1564.         return (intersects);
  1565.     else
  1566.         return behind;
  1567.     } else {
  1568.     if ((max_sign == 0) && (min_sign == 0))
  1569.         return (inside);
  1570.     else
  1571.         return infront;
  1572.     }
  1573. }
  1574.  
  1575.  
  1576. /* full xy-overlap test (support routines first) */
  1577.  
  1578. /* What do negative return values mean?
  1579.    It is allowed, that the 2 polygons touch themselves in one point.
  1580.    There are 2 possibilities of this case:
  1581.    (A) A vertex of the one polygon is on an edge of the other
  1582.    (B) The two polygons have a vertex together.
  1583.    In case (A) the algorithm finds 2 edges, wich touches an edge of the
  1584.    other polygon. In case (B) the algorithm finds four pairs of edges.
  1585.    I handle both cases with negative return values:
  1586.    Case (A) returns -2, and case (B) returns -1.
  1587.    So I can say, if the sum of negative return values goes greater than 4
  1588.    (absolutly), there must be more than one touch point. So the 
  1589.    polygons must overlap. */
  1590.  
  1591. /* some variables, for keeping track of minima and maxima across
  1592.  * all these routines */
  1593.  
  1594. static double m1x, M1x, m1y, M1y, m2x, M2x, m2y, M2y;
  1595.  
  1596. #define MINIMAX  (GE(M2x,m1x) && GE(M1x,m2x) && GE(M2y,m1y) && GE(M1y,m2y))
  1597.  
  1598.  
  1599. /* Does the edge [v1,v2] intersect edge [v3, v4] ?
  1600.  * This is for non-vertical [v1,v2]
  1601.  */
  1602. static int intersect(v1, v2, v3, v4)
  1603. struct Vertex GPHUGE *v1, GPHUGE * v2, GPHUGE * v3, GPHUGE * v4;
  1604. {
  1605.     double m1, m2, t1, t2, x, y, minx, maxx;
  1606.  
  1607.     m1 = (v2->y - v1->y) / (v2->x - v1->x);
  1608.     t1 = v1->y - m1 * v1->x;
  1609.  
  1610.     /* if [v3,v4] vertical: */
  1611.     if (EQ(v3->x, v4->x)) {
  1612.     y = v3->x * m1 + t1;
  1613.     if (GR(m2x, m1x) && GR(M1x, m2x)) {
  1614.         if (GR(y, m2y) && GR(M2y, y))
  1615.         return 1;
  1616.         if (EQ(y, m2y) || EQ(y, M2y))
  1617.         return -2;
  1618.         return 0;
  1619.     }
  1620.     /* m2x == m1x || m2x == M1x */
  1621.     if (GR(y, m2y) && GR(M2y, y))
  1622.         return -2;
  1623.     if (EQ(y, m2y) || EQ(y, M2y))
  1624.         return -1;
  1625.     return 0;
  1626.     }
  1627.     /* [v3,v4] not vertical */
  1628.     m2 = (v4->y - v3->y) / (v4->x - v3->x);
  1629.     t2 = v3->y - m2 * v3->x;
  1630.     if (!SIGNOF(m1 - m2)) {    /* if edges are parallel: */
  1631.     x = m1 * v3->x + t1 - v3->y;
  1632.     if (!EQ(m1 * v3->x + t1 - v3->y, 0.0))
  1633.         return 0;
  1634.     if (GR(M2x, m1x) && GR(M1x, m2x))
  1635.         return 1;
  1636.     return -1;        /* the edges have a common vertex */
  1637.     }
  1638.     x = (t2 - t1) / (m1 - m2);
  1639.     minx = GPMAX(m1x, m2x);
  1640.     maxx = GPMIN(M1x, M2x);
  1641.     if (GR(x, minx) && GR(maxx, x))
  1642.     return 1;
  1643.     if (GR(minx, x) || GR(x, maxx))
  1644.     return 0;
  1645.     if ((EQ(x, m1x) || EQ(x, M1x)) && (EQ(x, m1x) || EQ(x, M1x)))
  1646.     return -1;
  1647.     return -2;
  1648. }
  1649.  
  1650. /* Does the edge [v1,v2] intersect edge [v3, v4] ?
  1651.  * This is for vertical [v1,v2]
  1652.  */
  1653. static int v_intersect(v1, v2, v3, v4)
  1654. struct Vertex GPHUGE *v1, GPHUGE * v2, GPHUGE * v3, GPHUGE * v4;
  1655. {
  1656.     double y;
  1657.  
  1658.     /* HBB 971115: to silence picky compilers, use v2 at least once, in
  1659.      * a dummy expression (caution, may be ANSI-C specific because of
  1660.      * the use of 'void'...) */
  1661.     (void) v2;
  1662.  
  1663.     /* if [v3,v4] is vertical, too: */
  1664.     /* already known: rectangles do overlap, because MINIMAX was true */
  1665.     if (EQ(v3->x, v4->x))
  1666.     return (GR(M2y, m1y) && GR(M1y, m2y)) ? 1 : -1;
  1667.  
  1668.     /*  [v3,v4] not vertical */
  1669.     y = v3->y + (v1->x - v3->x) * (v4->y - v3->y) / (v4->x - v3->x);
  1670.     if (GR(m1x, m2x) && GR(M2x, m1x)) {
  1671.     if (GR(y, m1y) && GR(M1y, y))
  1672.         return 1;
  1673.     if (EQ(y, m1y) || EQ(y, M1y))
  1674.         return -2;
  1675.     return 0;
  1676.     }
  1677.     /* m1x == m2x || m1x == M2x */
  1678.     if (GR(y, m1y) && GR(M1y, y))
  1679.     return -2;
  1680.     if (EQ(y, m1y) || EQ(y, M1y))
  1681.     return -1;
  1682.     return 0;
  1683. }
  1684.  
  1685. #define UPD_MINMAX(v1,v2,npar) do {    \
  1686.   if (v1->x< v2->x)                 \
  1687.     CONCAT3(m,npar,x) = v1->x, CONCAT3(M,npar,x) = v2->x;   \
  1688.   else                              \
  1689.     CONCAT3(m,npar,x) = v2->x, CONCAT3(M,npar,x) = v1->x;   \
  1690.   if (v1->y< v2->y)                 \
  1691.     CONCAT3(m,npar,y) = v1->y, CONCAT3(M,npar,y) = v2->y;   \
  1692.   else                              \
  1693.     CONCAT3(m,npar,y) = v2->y, CONCAT3(M,npar,y) = v1->y;   \
  1694. } while (0)
  1695.  
  1696. /* does the edge [v1,v2] intersect polygon p? */
  1697. static int intersect_polygon(v1, v2, p)
  1698. struct Vertex GPHUGE *v1, GPHUGE * v2;
  1699. struct Polygon GPHUGE *p;
  1700. {
  1701.     struct Vertex GPHUGE *v3, GPHUGE * v4;
  1702.     int i, s, t = 0;
  1703.     int (*which_intersect) __PROTO((struct Vertex GPHUGE * v1,
  1704.                     struct Vertex GPHUGE * v2,
  1705.                     struct Vertex GPHUGE * v3,
  1706.                     struct Vertex GPHUGE * v4))
  1707.     = intersect;
  1708.  
  1709.     /* Is [v1,v2] vertical? If, use v_intersect() */
  1710.     if (EQ(v1->x, v2->x))
  1711.     which_intersect = v_intersect;
  1712.  
  1713.     UPD_MINMAX(v1,v2,1);
  1714.  
  1715.     /* test first edge of polygon p */
  1716.     v3 = vlist + p->vertex[p->n - 1];
  1717.     v4 = vlist + p->vertex[0];
  1718.     UPD_MINMAX(v3,v4,2);
  1719.  
  1720.  
  1721.     if (MINIMAX) {
  1722.     s = (*which_intersect) (v1, v2, v3, v4);
  1723.     if (s == 1 || (s < 0 && (t += s) < -4))
  1724.         return 1;
  1725.     }
  1726.     /* and the other edges... */
  1727.     for (i = 0; i < p->n - 1; i++) {
  1728.     v3 = vlist + p->vertex[i];
  1729.     v4 = vlist + p->vertex[i + 1];
  1730.     UPD_MINMAX(v3,v4,2);
  1731.     if (!MINIMAX)
  1732.         continue;
  1733.     s = (*which_intersect) (v1, v2, v3, v4);
  1734.     if (s == 1 || (s < 0 && (t += s) < -4))
  1735.         return 1;
  1736.     }
  1737.     return t;
  1738. }
  1739.  
  1740. /* OK, now here comes the 'real' polygon intersection routine:
  1741.  * Do a and b overlap in x or y (full test):
  1742.  * Do edges intersect, do the polygons touch in two points,
  1743.  * or is a vertex of one polygon inside the other polygon? */
  1744.  
  1745. static int full_xy_overlap(a, b)
  1746. struct Polygon GPHUGE *a, GPHUGE * b;
  1747. {
  1748.     struct Vertex GPHUGE *v1, GPHUGE * v2, v;
  1749.     int i, s, t = 0;
  1750.  
  1751.     if (a->n < 2 || b->n < 2)
  1752.     return 0;
  1753.     v1 = vlist + a->vertex[a->n - 1];
  1754.     v2 = vlist + a->vertex[0];
  1755.     s = intersect_polygon(v1, v2, b);
  1756.     if (s == 1 || (s < 0 && (t += s) < -4))
  1757.     return 1;
  1758.     for (i = 0; i < a->n - 1; i++) {
  1759.     v1 = vlist + a->vertex[i];
  1760.     v2 = vlist + a->vertex[i + 1];
  1761.     s = intersect_polygon(v1, v2, b);
  1762.     if (s == 1 || (s < 0 && (t += s) < -4))
  1763.         return 1;
  1764.     }
  1765.     /* No edges intersect. Is one polygon inside the other? */
  1766.     /* The  inner polygon has the greater minx */
  1767.     if (a->xmin < b->xmin) {
  1768.     struct Polygon GPHUGE *temp = a;
  1769.     a = b;
  1770.     b = temp;
  1771.     }
  1772.     /* Now, a is the inner polygon */
  1773.     for (i = 0; i < a->n; i++) {
  1774.     v1 = vlist + a->vertex[i];
  1775.     /* HBB: lclint seems to have found a problem here: v wasn't
  1776.      * fully defined when passed to intersect_polygon() */
  1777.     v = *v1;
  1778.     v.x = -1.1;
  1779.     s = intersect_polygon(v1, &v, b);
  1780.     if (s == 0)
  1781.         return 0;        /* a is outside polygon b */
  1782.     v.x = 1.1;
  1783.     if (s == 1) {
  1784.         s = intersect_polygon(v1, &v, b);
  1785.         if (s == 0)
  1786.         return 0;
  1787.         if (s == 1)
  1788.         return 1;
  1789.     } else if (intersect_polygon(v1, &v, b) == 0)
  1790.         return 0;
  1791.     }
  1792. #if 1
  1793.     /* 970614: don't bother, just presume that they do overlap: */
  1794.     return 1;
  1795. #else /* !1 */
  1796.     print_polygon(a, "a");
  1797.     print_polygon(b, "b");
  1798.     graph_error("Logic Error in full_xy_overlap");
  1799.     return -1;            /* HBB: shut up gcc -Wall */
  1800. #endif /* !1 */
  1801. }
  1802.  
  1803.  
  1804. /* Helper routine for split_polygon_by_plane */
  1805. static long build_new_vertex(V1, V2, w)
  1806. long V1, V2;    /* the vertices, between which a new vertex is demanded */
  1807. double w;    /* information about where between V1 and V2 it should be */
  1808. {
  1809.     long V;
  1810.     struct Vertex GPHUGE *v, GPHUGE * v1, GPHUGE * v2;
  1811.  
  1812.     if (EQ(w, 0.0))
  1813.     return V1;
  1814.     if (EQ(w, 1.0))
  1815.     return V2;
  1816.  
  1817.     /* We need a new Vertex */
  1818.     if (vert_free == nvert)    /* Extend vlist, if necessary */
  1819.     vlist = (struct Vertex GPHUGE *)
  1820.         gp_realloc(vlist, (unsigned long) sizeof(struct Vertex) * (nvert += 10L), "hidden vlist");
  1821.     V = vert_free++;
  1822.     v = vlist + V;
  1823.     v1 = vlist + V1;
  1824.     v2 = vlist + V2;
  1825.  
  1826.     v->x = (v2->x - v1->x) * w + v1->x;
  1827.     v->y = (v2->y - v1->y) * w + v1->y;
  1828.     v->z = (v2->z - v1->z) * w + v1->z;
  1829.     v->style = -1;
  1830.  
  1831.     /* 970610: additional checks, to prevent adding unnecessary vertices */
  1832.     if (V_EQUAL(v, v1)) {
  1833.     vert_free--;
  1834.     return V1;
  1835.     }
  1836.     if (V_EQUAL(v, v2)) {
  1837.     vert_free--;
  1838.     return V2;
  1839.     }
  1840.     return V;
  1841. }
  1842.  
  1843.  
  1844. /* Splits polygon p by the plane represented by its equation
  1845.  * coeffecients a to d.
  1846.  * return-value: part of the polygon, that is in front/back of plane
  1847.  * (Distinction necessary to ensure the ordering of polygons in
  1848.  * the plist stays intact)
  1849.  * Caution: plist and vlist can change their location!!!
  1850.  * If a point is in plane, it comes to the behind - part
  1851.  * (HBB: that may well have changed, as I cut at the 'in plane'
  1852.  * mechanisms heavily)
  1853.  */
  1854.  
  1855. static long split_polygon_by_plane(P, a, b, c, d, f)
  1856. long P;                /* Polygon as index on plist */
  1857. double a, b, c, d;
  1858. TBOOLEAN f;            /* return value = Front(1) or Back(0) */
  1859. {
  1860.     int i, j;
  1861.     struct Polygon GPHUGE *p = plist + P;
  1862.     struct Vertex GPHUGE *v;
  1863.     int snew, stest;
  1864.     int in_plane;        /* what is about points in plane? */
  1865.     int cross1, cross2;        /* first vertices, after p crossed the plane */
  1866.     double a1 = 0.0, b1, a2 = 0.0, b2;    /* parameters of the crosses */
  1867.     long Front;            /* the new Polygon */
  1868.     int n1, n2;
  1869.     long *v1, *v2;
  1870.     long line1, line2;
  1871.     long next_frag_temp;
  1872.  
  1873.  
  1874.     CHECK_POINTER(plist, p);
  1875.  
  1876.     in_plane = (EQ(c, 0.0) && f) ? 1 : -1;
  1877.     /* first vertex */
  1878.     v = vlist + p->vertex[0];
  1879.  
  1880.     CHECK_POINTER(vlist, v);
  1881.  
  1882.     b1 = a * v->x + b * v->y + c * v->z + d;
  1883.     stest = SIGNOF(b1);
  1884. #if OLD_SPLIT_INPLANE
  1885.     if (!stest)
  1886.     stest = in_plane;
  1887. #endif
  1888.  
  1889.     /* find first vertex that is on the other side of the plane */
  1890.     for (cross1 = 1, snew = stest; snew == stest && cross1 < p->n; cross1++) {
  1891.     a1 = b1;
  1892.     v = vlist + p->vertex[cross1];
  1893.     CHECK_POINTER(vlist, v);
  1894.     b1 = a * v->x + b * v->y + c * v->z + d;
  1895.     snew = SIGNOF(b1);
  1896. #if OLD_SPLIT_INPLANE
  1897.     if (!snew)
  1898.         snew = in_plane;
  1899. #endif
  1900.     }
  1901.  
  1902.     if (snew == stest) {
  1903.     /*HBB: this is possibly not an error in split_polygon, it's just
  1904.      * not possible to split this polygon by this plane. I.e., the
  1905.      * error is in the caller!  */
  1906.     fprintf(stderr, "split_poly failed, polygon nr. %ld\n", P);
  1907.     return (-1);        /* return 'unsplittable' */
  1908.     }
  1909.     cross1--;            /* now, it's the last one on 'this' side */
  1910.  
  1911.     /* first vertex that is on the first side again */
  1912.     for (b2 = b1, cross2 = cross1 + 1;
  1913.      snew != stest && cross2 < p->n; cross2++) {
  1914.     a2 = b2;
  1915.     v = vlist + p->vertex[cross2];
  1916.     CHECK_POINTER(vlist, v);
  1917.     b2 = a * v->x + b * v->y + c * v->z + d;
  1918.     snew = SIGNOF(b2);
  1919. #if OLD_SPLIT_INPLANE
  1920.     if (!snew)
  1921.         snew = -in_plane;
  1922. #endif
  1923.     }
  1924.  
  1925.     if (snew != stest) {
  1926.     /* only p->vertex[0] is on 'this' side */
  1927.     a2 = b2;
  1928.     v = vlist + p->vertex[0];
  1929.     CHECK_POINTER(vlist, v);
  1930.     b2 = a * v->x + b * v->y + c * v->z + d;
  1931.     } else
  1932.     cross2--;        /* now it's the last one on 'the other' side */
  1933.  
  1934.     /* We need two new polygons instead of the old one: */
  1935.     n1 = p->n - cross2 + cross1 + 2;
  1936.     n2 = cross2 - cross1 + 2;
  1937.     v1 = (long *) gp_alloc((unsigned long) sizeof(long) * n1,
  1938.                "hidden v1 for two new poly");
  1939.     v2 = (long *) gp_alloc((unsigned long) sizeof(long) * n2,
  1940.                "hidden v2 for two new poly");
  1941. #if SHOW_SPLITTING_EDGES
  1942.     line1 = 1L << (n1 - 1);
  1943.     line2 = 1L << (n2 - 1);
  1944. #else
  1945.     line1 = line2 = 0L;
  1946. #endif
  1947.     v1[0] = v2[n2 - 1] =
  1948.     build_new_vertex(p->vertex[cross2 - 1],
  1949.         p->vertex[(cross2 < p->n) ? cross2 : 0], a2 / (a2 - b2));
  1950.     v2[0] = v1[n1 - 1] =
  1951.     build_new_vertex(p->vertex[cross1 - 1],
  1952.              p->vertex[cross1], a1 / (a1 - b1));
  1953.  
  1954.     /* Copy visibility from split edges to their fragments: */
  1955.     if (p->line & (1 << (cross1 - 1))) {
  1956.     line1 |= 1L << (n1 - 2);
  1957.     line2 |= 1L;
  1958.     }
  1959.     if (p->line & (1L << (cross2 - 1))) {
  1960.     line1 |= 1L;
  1961.     line2 |= 1L << (n2 - 2);
  1962.     }
  1963.     /* Copy the old line patterns and vertex listings into new places */
  1964.     /* p->vertex[cross2...p->n-1] --> v1[1...] */
  1965.     for (i = cross2, j = 1; i < p->n;) {
  1966.     if (p->line & (1L << i))
  1967.         line1 |= 1L << j;
  1968.     v1[j++] = p->vertex[i++];
  1969.     }
  1970.     /* p->vertex[0...cross1-1] --> v1[...n1-2] */
  1971.     for (i = 0; i < cross1 - 1;) {
  1972.     if (p->line & (1L << i))
  1973.         line1 |= 1L << j;
  1974.     v1[j++] = p->vertex[i++];
  1975.     }
  1976.     v1[j++] = p->vertex[i++];
  1977.  
  1978.     /* p->vertex[cross1...cross2-1] -> v2[1...n2-2] */
  1979.     for (j = 1; i < cross2 - 1;) {
  1980.     if (p->line & (1L << i))
  1981.         line2 |= 1L << j;
  1982.     v2[j++] = p->vertex[i++];
  1983.     }
  1984.     v2[j++] = p->vertex[i];
  1985.  
  1986.     /* old vertex list isn't needed any more: */
  1987.     free(p->vertex);
  1988.     p->vertex = 0;
  1989.  
  1990.     if (!reduce_polygon(&n1, &v1, &line1, 1) ||
  1991.     !reduce_polygon(&n2, &v2, &line2, 1))
  1992.     graph_error("Logic Error 2 in split_polygon");
  1993.  
  1994.     /* Build the 2 new polygons : we are reusing one + making one new */
  1995.     CHECK_PLIST_EXTRA(p = plist + P);
  1996.     Front = pfree++;
  1997.  
  1998.     if ((next_frag_temp = p->next_frag) < 0)
  1999.     next_frag_temp = P;    /* First split of this polygon at all */
  2000.  
  2001.     /* HBB 961110: lclint said == shouldn't be used on Boolean's */
  2002.     if ((f && (stest < 0)) || ((!f) && !(stest < 0))) {
  2003.     build_polygon(plist + P, n1, v1, line1, p->style, p->lp,
  2004.               p->next, Front, p->id, p->tested);
  2005.     build_polygon(plist + Front, n2, v2, line2, p->style, p->lp,
  2006.               -1, next_frag_temp, p->id, is_moved_or_split);
  2007.     } else {
  2008.     build_polygon(plist + P, n2, v2, line2, p->style, p->lp,
  2009.               p->next, Front, p->id, p->tested);
  2010.     build_polygon(plist + Front, n1, v1, line1, p->style, p->lp,
  2011.               -1, next_frag_temp, p->id, is_moved_or_split);
  2012.     }
  2013.     return Front;
  2014. }
  2015.  
  2016. /* HBB: these pieces of code are found repeatedly in in_front(), so
  2017.  * I put them into macros
  2018.  * Note: can't use the 'do {st, This))
  2019.         method for
  2020.  * these: 'continue' wouldn't work any more.
  2021.  */
  2022. #define PUT_P_IN_FRONT_TEST(new_tested) {\
  2023.   plist[Plast].next = p->next; \
  2024.   p->next = Test; \
  2025.   p->tested = new_tested; \
  2026.   if (Insert >= 0) \
  2027.     plist[Insert].next = P; \
  2028.   else \
  2029.     pfirst = P; \
  2030.   Insert = P; \
  2031.   P = Plast; \
  2032.   p = plist + P; \
  2033.   continue; \
  2034. }
  2035.  
  2036. #define SPLIT_TEST_BY_P {\
  2037.   long Back = split_polygon_by_plane (Test, pa, pb, pc, pd, 0); \
  2038.   if (Back >0) { \
  2039.     p = plist + P; \
  2040.     test = plist + Test;  /* plist can change */ \
  2041.     plist[Back].next = p->next; \
  2042.     p->next = Back; \
  2043.     P = Back; \
  2044.     p = plist + P; \
  2045.     zmin = test->zmin;  /* the new z-value */ \
  2046.   } else { \
  2047.     fprintf(stderr, "Failed to split test (%ld) by p (%ld)\n", Test, P); \
  2048.     graph_error("Error in hidden line removal: couldn't split..."); \
  2049.   } \
  2050.   continue; \
  2051. }
  2052.  
  2053. #define SPLIT_P_BY_TEST {\
  2054.   long Front = split_polygon_by_plane (P, ta, tb, tc, td, 1);\
  2055.   if (Front >0) {\
  2056.     p = plist + P;\
  2057.     test = plist + Test;    /* plist can change */\
  2058.     plist[Front].next = Test;\
  2059.     if (Insert >= 0)\
  2060.       plist[Insert].next = Front;\
  2061.     else\
  2062.       pfirst = Front;\
  2063.     Insert = Front;\
  2064.   } else {\
  2065.     fprintf(stderr, "Failed to split p (%ld) by test(%ld), relations are %d, %d\n", P, Test, p_rel_tplane, t_rel_pplane);\
  2066.     print_polygon(test, "test");\
  2067.     print_polygon(p, "p");\
  2068.     fputc('\n', stderr);\
  2069.   }\
  2070.   continue; /* FIXME: should we continue nevertheless? */\
  2071. }
  2072.  
  2073.  
  2074. /* Is Test in front of the other polygons?
  2075.  * If not, polygons which are in front of test come between
  2076.  * Last and Test (I.e. to the 'front' of the plist)
  2077.  */
  2078. static int in_front(Last, Test)
  2079. long Last, Test;
  2080. {
  2081.     struct Polygon GPHUGE *test = plist + Test, GPHUGE * p;
  2082.     long Insert, P, Plast;
  2083.     double zmin,        /* minimum z-value of Test */
  2084.      ta, tb, tc, td,        /* the plane of Test      */
  2085.      pa, pb, pc, pd;        /* the plane of a polygon */
  2086.     t_poly_plane_intersect p_rel_tplane, t_rel_pplane;
  2087.     int loop = (test->tested == is_in_loop);
  2088.  
  2089.     CHECK_POINTER(plist, test);
  2090.  
  2091.     /* Maybe this should only done at the end of this routine... */
  2092.     /* test->tested = is_tested; */
  2093.  
  2094.     Insert = Last;
  2095.  
  2096.     /* minimal z-value of Test */
  2097.     zmin = test->zmin;
  2098.  
  2099.     /* The plane of the polygon Test */
  2100.     get_plane(test, &ta, &tb, &tc, &td);
  2101.  
  2102.     /* Compare Test with the following polygons, which overlap in z value */
  2103.     for (Plast = Test, p = test;
  2104.      ((P = p->next) >= 0) && (((p = plist + P)->zmax > zmin) || (p->tested != 0));
  2105.      Plast = P) {
  2106.     CHECK_POINTER(plist, p);
  2107.  
  2108.     if (!xy_overlap(test, p))
  2109.         continue;
  2110.  
  2111.     if (p->zmax <= zmin)
  2112.         continue;
  2113.  
  2114.     p_rel_tplane = polygon_plane_intersection(p, ta, tb, tc, td);
  2115.     if ((p_rel_tplane == behind) || (p_rel_tplane == inside))
  2116.         continue;
  2117.  
  2118.     get_plane(p, &pa, &pb, &pc, &pd);
  2119.     t_rel_pplane = polygon_plane_intersection(test, pa, pb, pc, pd);
  2120.     if ((t_rel_pplane == infront) || (t_rel_pplane == inside))
  2121.         continue;
  2122.  
  2123.     if (!full_xy_overlap(test, p))
  2124.         continue;
  2125.  
  2126. #if 1
  2127.     /* HBB 971115: try new approach developed directly from Foley et
  2128.      * al.'s original description */
  2129.     /* OK, at this point, a total of 16 cases is left to be handled,
  2130.      * as 4 variables are important, each with two possible values:
  2131.      * t_rel_pplane may be 'behind' or 'intersects', p_rel_tplane may
  2132.      * be 'infront' or 'intersects', 'test' may be (seem to be) part
  2133.      * of a loop or not ('loop'), and p may have been considered
  2134.      * already or not (p->tested =?= is_tested). */
  2135.     /* See if splitting wherever possible is a good idea: */
  2136.     if (p_rel_tplane == intersects) {
  2137.         p->tested = is_moved_or_split;
  2138.         SPLIT_P_BY_TEST;
  2139.     } else if (t_rel_pplane == intersects) {
  2140.         test->tested = is_moved_or_split;
  2141.         SPLIT_TEST_BY_P;
  2142.     } else {
  2143.         if (loop && (p->tested == is_tested)) {
  2144.         /* Ouch, seems like we're in trouble, really */
  2145.         fprintf(stderr, "\
  2146. #Failed: In loop, without intersections.\n\
  2147. #Relations are %d, %d\n",
  2148.             p_rel_tplane, t_rel_pplane);
  2149.         print_polygon(test, "test");
  2150.         print_polygon(p, "p");
  2151.         continue;    /* Keep quiet, maybe no-one will notice (;-) */
  2152.         } else {
  2153.         PUT_P_IN_FRONT_TEST(is_in_loop);
  2154.         }
  2155.     }
  2156. #else
  2157.     /* p obscures test (at least, this *might* be happening */
  2158.     if (loop) {
  2159.         /* if P isn't interesting for the loop, put P in front of Test */
  2160.         if ((p->tested != is_tested)
  2161.         || (p_rel_tplane == infront)
  2162.         || (t_rel_pplane == behind))    /* HBB 970325: was infront (!?) */
  2163.         PUT_P_IN_FRONT_TEST(is_moved_or_split);
  2164.  
  2165.         /* else, split either test or p */
  2166.         if (t_rel_pplane == intersects)
  2167.         SPLIT_TEST_BY_P;
  2168.         if (p_rel_tplane == intersects)
  2169.         SPLIT_P_BY_TEST;
  2170.  
  2171.         /* HBB: if we ever come through here, something went severely wrong */
  2172.         fprintf(stderr, "Failed: polygon %ld vs. %ld, relations are %d, %d\n", P, Test, p_rel_tplane, t_rel_pplane);
  2173.         print_polygon(test, "test");
  2174.         print_polygon(p, "p");
  2175.         fputc('\n', stderr);
  2176.         graph_error("Couldn't resolve a polygon overlapping pb. Go tell HBB...");
  2177.     }
  2178.     /* No loop: if it makes sense, put P in front of Test */
  2179.     if ((p->tested != is_tested)
  2180.         && ((t_rel_pplane == behind)
  2181.         || (p_rel_tplane == infront)))
  2182.         PUT_P_IN_FRONT_TEST(is_moved_or_split);
  2183.  
  2184.     /* If possible, split P */
  2185.     if ((p->tested != is_tested) || (p_rel_tplane == intersects))
  2186.         SPLIT_P_BY_TEST;
  2187.  
  2188.     /* if possible, split Test */
  2189.     if (t_rel_pplane == intersects)
  2190.         SPLIT_TEST_BY_P;
  2191.  
  2192.     /* else, we have a loop: mark P as loop and put it in front of Test */
  2193.     PUT_P_IN_FRONT_TEST(is_in_loop);    /* -2: p is in a loop */
  2194. #endif
  2195.     }
  2196.     if (Insert == Last)
  2197.     test->tested = is_tested;
  2198.     return (int) (Insert == Last);
  2199. }
  2200.  
  2201.  
  2202. /* Drawing the polygons */
  2203.  
  2204. /* HBB 980303: new functions to handle Cross back-buffer (saves
  2205.  * gp_alloc()'s) */
  2206.  
  2207. static GP_INLINE void init_Cross_store()
  2208. {
  2209.     assert(!Cross_store && !last_Cross_store);
  2210.     last_Cross_store = CROSS_STORE_INCREASE;
  2211.     free_Cross_store = 0;
  2212.     Cross_store = (struct Cross *)
  2213.     gp_alloc((unsigned long) last_Cross_store * sizeof(struct Cross),
  2214.          "hidden cross store");
  2215. }
  2216.  
  2217. static GP_INLINE struct Cross *get_Cross_from_store()
  2218. {
  2219.     while (last_Cross_store <= free_Cross_store) {
  2220.     last_Cross_store += CROSS_STORE_INCREASE;
  2221.     Cross_store = (struct Cross *)
  2222.         gp_realloc(Cross_store,
  2223.          (unsigned long) last_Cross_store * sizeof(struct Cross),
  2224.                "hidden cross store");
  2225.     }
  2226.     return Cross_store + (free_Cross_store++);
  2227. }
  2228.  
  2229. static GP_INLINE TBOOLEAN
  2230.  hl_buffer_set(xv, yv)
  2231. int xv, yv;
  2232. {
  2233.     struct Cross GPHUGE *c;
  2234.     /*HBB 961110: lclint wanted this: */
  2235.     assert(hl_buffer != 0);
  2236.     for (c = hl_buffer[xv]; c != NULL; c = c->next)
  2237.     if (c->a <= yv && c->b >= yv) {
  2238.         return TRUE;
  2239.     }
  2240.     return FALSE;
  2241. }
  2242.  
  2243. /* HBB 961201: new var's, to speed up free()ing hl_buffer later */
  2244. static int hl_buff_xmin, hl_buff_xmax;
  2245.  
  2246. /* HBB 961124: new routine. All this occured twice in
  2247.  * draw_clip_line_buffer */
  2248. /* Store a line crossing the x interval around xv between y = ya and
  2249.  * y = yb in the hl_buffer */
  2250. static GP_INLINE void update_hl_buffer_column(xv, ya, yb)
  2251. int xv, ya, yb;
  2252. {
  2253.     struct Cross GPHUGE *GPHUGE * cross, GPHUGE * cross2;
  2254.  
  2255.     /* First, ensure that ya <= yb */
  2256.     if (ya > yb) {
  2257.     int y_temp = yb;
  2258.     yb = ya;
  2259.     ya = y_temp;
  2260.     }
  2261.     /* loop over all previous crossings at this x-value */
  2262.     for (cross = hl_buffer + xv; 1; cross = &(*cross)->next) {
  2263.     if (*cross == NULL) {
  2264.         /* first or new highest crossing at this x-value */
  2265.  
  2266.         /* HBB 980303: new method to allocate Cross structures */
  2267.         *cross = get_Cross_from_store();
  2268.         (*cross)->a = ya;
  2269.         (*cross)->b = yb;
  2270.         (*cross)->next = NULL;
  2271.         /* HBB 961201: keep track of x-range of hl_buffer, to
  2272.          * speedup free()ing it */
  2273.         if (xv < hl_buff_xmin)
  2274.         hl_buff_xmin = xv;
  2275.         if (xv > hl_buff_xmax)
  2276.         hl_buff_xmax = xv;
  2277.         break;
  2278.     }
  2279.     if (yb < (*cross)->a - 1) {
  2280.         /* crossing below 'cross', create new entry before 'cross' */
  2281.         cross2 = *cross;
  2282.         /* HBB 980303: new method to allocate Cross structures */
  2283.         *cross = get_Cross_from_store();
  2284.         (*cross)->a = ya;
  2285.         (*cross)->b = yb;
  2286.         (*cross)->next = cross2;
  2287.         break;
  2288.     } else if (ya <= (*cross)->b + 1) {
  2289.         /* crossing overlaps or covers 'cross' */
  2290.         if (ya < (*cross)->a)
  2291.         (*cross)->a = ya;
  2292.         if (yb > (*cross)->b) {
  2293.         if ((*cross)->next && (*cross)->next->a <= yb) {
  2294.             /* crossing spans all the way up to 'cross->next' so
  2295.              * unite them */
  2296.             cross2 = (*cross)->next;
  2297.             (*cross)->b = cross2->b;
  2298.             (*cross)->next = cross2->next;
  2299.         } else
  2300.             (*cross)->b = yb;
  2301.         }
  2302.         break;
  2303.     }
  2304.     }
  2305.     return;
  2306. }
  2307.  
  2308.  
  2309. static void draw_clip_line_buffer(x1, y1, x2, y2)
  2310. int x1, y1, x2, y2;
  2311.      /* Draw a line in the hl_buffer */
  2312. {
  2313.     register int xv, yv, errx, erry, err;
  2314.     register int dy, nstep;
  2315.     register int ya;
  2316.     int i;
  2317.  
  2318.     if (!clip_line(&x1, &y1, &x2, &y2)) {
  2319.     /*printf("dcl_buffer: clipped away!\n"); */
  2320.     return;
  2321.     }
  2322.     if (x1 > x2) {
  2323.     errx = XREDUCE(x1) - XREDUCE(x2);
  2324.     erry = YREDUCE(y1) - YREDUCE(y2);
  2325.     xv = XREDUCE(x2) - XREDUCE(xleft);
  2326.     yv = YREDUCE(y2) - YREDUCE(ybot);
  2327.     } else {
  2328.     errx = XREDUCE(x2) - XREDUCE(x1);
  2329.     erry = YREDUCE(y2) - YREDUCE(y1);
  2330.     xv = XREDUCE(x1) - XREDUCE(xleft);
  2331.     yv = YREDUCE(y1) - YREDUCE(ybot);
  2332.     }
  2333.     if (erry > 0)
  2334.     dy = 1;
  2335.     else {
  2336.     dy = -1;
  2337.     erry = -erry;
  2338.     }
  2339.     nstep = errx + erry;
  2340.     err = -errx + erry;
  2341.     errx <<= 1;
  2342.     erry <<= 1;
  2343.     ya = yv;
  2344.  
  2345.     for (i = 0; i < nstep; i++) {
  2346.     if (err < 0) {
  2347.         update_hl_buffer_column(xv, ya, yv);
  2348.         xv++;
  2349.         ya = yv;
  2350.         err += erry;
  2351.     } else {
  2352.         yv += dy;
  2353.         err -= errx;
  2354.     }
  2355.     }
  2356.     (void) update_hl_buffer_column(xv, ya, yv);
  2357.     return;
  2358. }
  2359.  
  2360.  
  2361. /* HBB 961124: improve similarity of this routine with
  2362.  * draw_clip_line_buffer ()
  2363.  *
  2364.  * HBB 961124: changed from 'virtual' to 'do_draw', for clarity (this
  2365.  * also affects the routines calling this one, so beware!
  2366.  *
  2367.  * HBB 961124: introduced checking code, to search for the 'holes in
  2368.  * the surface' bug
  2369.  *
  2370.  * Draw a line into the bitmap (**cpnt), and optionally also to the
  2371.  * terminal
  2372.  */
  2373. static void draw_clip_line_update(x1, y1, x2, y2, do_draw)
  2374. int x1, y1, x2, y2;
  2375. TBOOLEAN do_draw;
  2376. {
  2377.     /* HBB 961110: made 'flag' a boolean variable, which needs some
  2378.      *  other changes below as well
  2379.      *
  2380.      * HBB 970508: renamed 'flag' to 'is_drawn' (reverting its meaning).
  2381.      * Should be clearer now
  2382.      */
  2383.     TBOOLEAN is_drawn;
  2384.     register struct termentry *t = term;
  2385.     register int xv, yv, errx, erry, err;
  2386.     register int xvr, yvr;
  2387.     int xve, yve;
  2388.     register int dy, nstep = 0, dyr;
  2389.     int i;
  2390.     if (!clip_line(&x1, &y1, &x2, &y2)) {
  2391.     /*printf("dcl_update: clipped away!\n"); */
  2392.     return;
  2393.     }
  2394.     if (x1 > x2) {
  2395.     xvr = x2;
  2396.     yvr = y2;
  2397.     xve = x1;
  2398.     yve = y1;
  2399.     } else {
  2400.     xvr = x1;
  2401.     yvr = y1;
  2402.     xve = x2;
  2403.     yve = y2;
  2404.     }
  2405.     errx = XREDUCE(xve) - XREDUCE(xvr);
  2406.     erry = YREDUCE(yve) - YREDUCE(yvr);
  2407.     if (erry > 0)
  2408.     dy = 1;
  2409.     else {
  2410.     dy = -1;
  2411.     erry = -erry;
  2412.     }
  2413.     dyr = dy * yfact;
  2414.     nstep = errx + erry;
  2415.     err = -errx + erry;
  2416.     errx <<= 1;
  2417.     erry <<= 1;
  2418.     xv = XREDUCE(xvr) - XREDUCE(xleft);
  2419.     yv = YREDUCE(yvr) - YREDUCE(ybot);
  2420.     (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2421.  
  2422.     is_drawn = (IS_UNSET(xv, yv) && (do_draw || hl_buffer_set(xv, yv)));
  2423.  
  2424.     /* HBB 961110: lclint wanted these: */
  2425.     assert(ymin_hl != 0);
  2426.     assert(ymax_hl != 0);
  2427.     /* Check first point in hl-buffer */
  2428.     if (xv < xmin_hl)
  2429.     xmin_hl = xv;
  2430.     if (xv > xmax_hl)
  2431.     xmax_hl = xv;
  2432.     if (yv < ymin_hl[xv])
  2433.     ymin_hl[xv] = yv;
  2434.     if (yv > ymax_hl[xv])
  2435.     ymax_hl[xv] = yv;
  2436.     for (i = 0; i < nstep; i++) {
  2437.     /* 970722: new variables */
  2438.     unsigned int xvr_temp = xvr, yvr_temp = yvr;
  2439.     if (err < 0) {
  2440.         xv++;
  2441.         xvr += xfact;
  2442.         err += erry;
  2443.     } else {
  2444.         yv += dy;
  2445.         yvr += dyr;
  2446.         err -= errx;
  2447.     }
  2448.     if (IS_UNSET(xv, yv) && (do_draw || hl_buffer_set(xv, yv))) {
  2449.         if (!is_drawn) {
  2450.         /* 970722: If this line was meant to be drawn, originally,
  2451.          * and it was just the bitmap that stopped it, then draw
  2452.          * it a bit longer (i.e.: start one pixel earlier
  2453.          */
  2454.         if (do_draw)
  2455.             (*t->move) ((unsigned int) xvr_temp, (unsigned int) yvr_temp);
  2456.         else
  2457.             (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2458.         is_drawn = TRUE;
  2459.         }
  2460.     } else {
  2461.         if (is_drawn) {
  2462.         /* 970722: If this line is only drawn because hl_buffer_set()
  2463.          * was true, then don't extend to the new position (where it
  2464.          * isn't true any more). Draw the line one pixel shorter,
  2465.          * instead:
  2466.          */
  2467.         if (do_draw)
  2468.             (*t->vector) ((unsigned int) xvr, (unsigned int) yvr);
  2469.         else {
  2470.             (*t->vector) ((unsigned int) xvr_temp, (unsigned int) yvr_temp);
  2471.             (*t->move) ((unsigned int) xvr, (unsigned int) yvr);
  2472.         }
  2473.         is_drawn = FALSE;
  2474.         }
  2475.     }
  2476.     if (xv < xmin_hl)
  2477.         xmin_hl = xv;
  2478.     if (xv > xmax_hl)
  2479.         xmax_hl = xv;
  2480.     if (yv < ymin_hl[xv])
  2481.         ymin_hl[xv] = yv;
  2482.     if (yv > ymax_hl[xv])
  2483.         ymax_hl[xv] = yv;
  2484.     }
  2485.     if (is_drawn)
  2486.     (*t->vector) ((unsigned int) xvr, (unsigned int) yvr);
  2487.     return;
  2488. }
  2489.  
  2490. #define DRAW_VERTEX(v, x, y) do { \
  2491.    if ((v)->style >= 0 &&                                            \
  2492.        !clip_point(x,y)  &&                                        \
  2493.        IS_UNSET(XREDUCE(x)-XREDUCE(xleft),YREDUCE(y)-YREDUCE(ybot))) \
  2494.      (*t->point)(x,y, v->style);                                   \
  2495.    (v)->style = -1;                                                  \
  2496. } while (0)
  2497.  
  2498. /* Two routines to emulate move/vector sequence using line drawing routine. */
  2499. static unsigned int move_pos_x, move_pos_y;
  2500.  
  2501. void clip_move(x, y)
  2502. unsigned int x, y;
  2503. {
  2504.     move_pos_x = x;
  2505.     move_pos_y = y;
  2506. }
  2507.  
  2508. void clip_vector(x, y)
  2509. unsigned int x, y;
  2510. {
  2511.     draw_clip_line(move_pos_x, move_pos_y, x, y);
  2512.     move_pos_x = x;
  2513.     move_pos_y = y;
  2514. }
  2515.  
  2516. static GP_INLINE void clip_vector_h(x, y)
  2517. int x, y;
  2518.      /* Draw a line on terminal and update bitmap */
  2519. {
  2520.     draw_clip_line_update(move_pos_x, move_pos_y, x, y, TRUE);
  2521.     move_pos_x = x;
  2522.     move_pos_y = y;
  2523. }
  2524.  
  2525.  
  2526. static GP_INLINE void clip_vector_virtual(x, y)
  2527. int x, y;
  2528.      /* update bitmap, do not really draw the line */
  2529. {
  2530.     draw_clip_line_update(move_pos_x, move_pos_y, x, y, FALSE);
  2531.     move_pos_x = x;
  2532.     move_pos_y = y;
  2533. }
  2534.  
  2535. static GP_INLINE void clip_vector_buffer(x, y)
  2536.      /* draw a line in the hl_buffer */
  2537. int x, y;
  2538. {
  2539.     draw_clip_line_buffer(move_pos_x, move_pos_y, x, y);
  2540.     move_pos_x = x;
  2541.     move_pos_y = y;
  2542. }
  2543.  
  2544. static void draw(p)
  2545. struct Polygon GPHUGE *p;
  2546. {
  2547.     struct Vertex GPHUGE *v;
  2548.     struct Polygon GPHUGE *q;
  2549.     long Q;
  2550.     int i;
  2551.     int x, y;            /* point in terminal coordinates */
  2552.     register struct termentry *t = term;
  2553.     t_pnt mask1, mask2;
  2554.     long int indx1, indx2, k;
  2555.     tp_pnt cpnt;
  2556.  
  2557.     xmin_hl = HL_EXTENT_X_MAX;    /* HBB 961124: parametrized this value */
  2558.     xmax_hl = 0;
  2559.  
  2560.     /* HBB 961201: store initial values for range of hl_buffer: */
  2561.     hl_buff_xmin = XREDUCE(xright) - XREDUCE(xleft);
  2562.     hl_buff_xmax = 0;
  2563.  
  2564.     /* Write the lines of all polygons with the same id as p in the hl_buffer */
  2565.     /* HBB 961201: use next_frag pointers to loop over fragments: */
  2566.     for (q = p; (Q = q->next_frag) >= 0 && (q = plist + Q) != p;) {
  2567.     if (q->id != p->id)
  2568.         continue;
  2569.  
  2570.     /* draw the lines of q into hl_buffer: */
  2571.     v = vlist + q->vertex[0];
  2572.     TERMCOORD(v, x, y);
  2573.     clip_move(x, y);
  2574.     for (i = 1; i < q->n; i++) {
  2575.         v = vlist + q->vertex[i];
  2576.         TERMCOORD(v, x, y);
  2577.         if (q->line & (1L << (i - 1)))
  2578.         clip_vector_buffer(x, y);
  2579.         else
  2580.         clip_move(x, y);
  2581.     }
  2582.     if (q->line & (1L << (i - 1))) {
  2583.         v = vlist + q->vertex[0];
  2584.         TERMCOORD(v, x, y);
  2585.         clip_vector_buffer(x, y);
  2586.     }
  2587.     }
  2588.  
  2589.  
  2590.     /* first, set the line and point styles: */
  2591.     {
  2592.     struct lp_style_type lp;
  2593.  
  2594.     lp = *(p->lp);
  2595.     lp.l_type = p->style;
  2596.     term_apply_lp_properties(&(lp));
  2597.     }
  2598.  
  2599.     /*print_polygon (p, "drawing"); */
  2600.  
  2601.     /* draw the lines of p */
  2602.     v = vlist + p->vertex[0];
  2603.     TERMCOORD(v, x, y);
  2604.     clip_move(x, y);
  2605.     DRAW_VERTEX(v, x, y);
  2606.     for (i = 1; i < p->n; i++) {
  2607.     v = vlist + p->vertex[i];
  2608.     TERMCOORD(v, x, y);
  2609.     if (p->line & (1L << (i - 1)))
  2610.         clip_vector_h(x, y);
  2611.     else
  2612.         clip_vector_virtual(x, y);
  2613.     DRAW_VERTEX(v, x, y);
  2614.     }
  2615.     TERMCOORD(vlist + p->vertex[0], x, y);
  2616.     if (p->line & (1L << (i - 1)))
  2617.     clip_vector_h(x, y);
  2618.     else
  2619.     clip_vector_virtual(x, y);
  2620.  
  2621.  
  2622.     /* reset the hl_buffer */
  2623.     /*HBB 961110: lclint wanted this: */
  2624.     assert(hl_buffer);
  2625.     for (i = hl_buff_xmin; i <= hl_buff_xmax; i++) {
  2626.     /* HBB 980303: one part was removed here. It isn't needed any more,
  2627.      * with the global store for Cross structs. */
  2628.     hl_buffer[i] = NULL;
  2629.     }
  2630.     /* HBB 980303: instead, set back the free pointer of the Cross store: */
  2631.     free_Cross_store = 0;
  2632.  
  2633.     /* now mark the area as being filled in the bitmap. */
  2634.     /* HBB 971115: Yes, I do know that xmin_hl is unsigned, and the
  2635.      * first test is thus useless. But I'd like to keep it anyway ;-) */
  2636.     if (xmin_hl < 0 || xmax_hl > XREDUCE(xright) - XREDUCE(xleft))
  2637.     graph_error("Logic error #3 in hidden line");
  2638.     /* HBB 961110: lclint wanted these: */
  2639.     assert(ymin_hl != 0);
  2640.     assert(ymax_hl != 0);
  2641.     assert(pnt != 0);
  2642.     if (xmin_hl < xmax_hl)
  2643.     for (i = xmin_hl; i <= xmax_hl; i++) {
  2644.         if (ymin_hl[i] == HL_EXTENT_Y_MAX)
  2645.         graph_error("Logic error #2 in hidden line");
  2646.         if (pnt[i] == 0) {
  2647.         pnt[i] = (t_pnt *) gp_alloc((unsigned long) y_malloc, "hidden ymalloc");
  2648.         memset(pnt[i], 0, (size_t) y_malloc);
  2649.         }
  2650.         if (ymin_hl[i] < 0 || ymax_hl[i] > YREDUCE(ytop) - YREDUCE(ybot))
  2651.         graph_error("Logic error #1 in hidden line");
  2652.         /* HBB 970508: this shift was wordsize dependent
  2653.          * I think I fixed that
  2654.          */
  2655.         indx1 = ymin_hl[i] / PNT_BITS;
  2656.         indx2 = ymax_hl[i] / PNT_BITS;
  2657.         mask1 = PNT_MAX << (((unsigned) ymin_hl[i]) % PNT_BITS);
  2658.         mask2 = PNT_MAX >> ((~((unsigned) ymax_hl[i])) % PNT_BITS);
  2659.         cpnt = pnt[i] + indx1;
  2660.         if (indx1 == indx2)
  2661.         *cpnt |= (mask1 & mask2);
  2662.         else {
  2663.         *cpnt++ |= mask1;
  2664.         for (k = indx1 + 1; k < indx2; k++)
  2665.             *cpnt++ = PNT_MAX;
  2666.         *cpnt |= mask2;
  2667.         }
  2668.         ymin_hl[i] = HL_EXTENT_Y_MAX;
  2669.         ymax_hl[i] = 0;
  2670.     }
  2671. }
  2672.  
  2673.  
  2674. void plot3d_hidden(plots, pcount)
  2675. struct surface_points *plots;
  2676. int pcount;
  2677. {
  2678.     long Last, This;
  2679.     long i;
  2680.  
  2681. #if defined(DOS16) || defined(WIN16)
  2682.     /* HBB 980309: Ensure that Polygon Structs exactly fit a 64K segment. The
  2683.      * problem this prevents is that even with 'huge pointers', a single
  2684.      * struct may not cross a 64K boundary: it will wrap around to the
  2685.      * beginning of that 64K segment. Someone ought to complain about
  2686.      * this nuisance, somewhere... :-( */
  2687.     assert(0 == (0x1 << 16) % (sizeof(struct Polygon)));
  2688. #endif
  2689.  
  2690.     /* Initialize global variables */
  2691.     y_malloc = (2 + (YREDUCE(ytop) >> 4) - (YREDUCE(ybot) >> 4)) * sizeof(t_pnt);
  2692.     /* ymin_hl, ymax_hl: */
  2693.     i = sizeof(t_hl_extent_y) * (XREDUCE(xright) - XREDUCE(xleft) + 1);
  2694.     ymin_hl = (t_hl_extent_y *) gp_alloc((unsigned long) i, "hidden ymin_hl");
  2695.     ymax_hl = (t_hl_extent_y *) gp_alloc((unsigned long) i, "hidden ymax_hl");
  2696.     for (i = (XREDUCE(xright) - XREDUCE(xleft)); i >= 0; i--) {
  2697.     ymin_hl[i] = HL_EXTENT_Y_MAX;
  2698.     ymax_hl[i] = 0;
  2699.     }
  2700.     /* hl_buffer: */
  2701.     /* HBB 980303 new: initialize the global store for Cross structs: */
  2702.     init_Cross_store();
  2703.     i = XREDUCE(xright) - XREDUCE(xleft) + 1;
  2704.     hl_buffer =
  2705.     (struct Cross GPHUGE * GPHUGE *) gp_alloc((unsigned long) (i * sizeof(struct Cross GPHUGE *)),
  2706.                           "hidden hl_buffer");
  2707.     while (--i >= 0)
  2708.     hl_buffer[i] = (struct Cross *) 0;
  2709.  
  2710.     init_polygons(plots, pcount);
  2711.  
  2712.     /* HBB 970326: if no polygons survived the cutting away of undefined
  2713.      * and/or outranged vertices, bail out with a warning message.
  2714.      * Without this, sort_by_zmax gets a SegFault */
  2715.     if (!pfree) {
  2716.     /* HBB 980701: new code to deallocate all the structures allocated up
  2717.      * to this point. Should fix the core dump reported by 
  2718.      * Stefan A. Deutscher today  */
  2719.     if (ymin_hl) {
  2720.         free(ymin_hl);
  2721.         ymin_hl = 0;
  2722.     }
  2723.     if (ymax_hl) {
  2724.         free(ymax_hl);
  2725.         ymax_hl = 0;
  2726.     }
  2727.     if (hl_buffer) {
  2728.         free(hl_buffer);
  2729.         hl_buffer = 0;
  2730.     }
  2731.     if (Cross_store) {
  2732.         free(Cross_store);
  2733.         Cross_store = 0;
  2734.         last_Cross_store = 0;
  2735.     }
  2736.     graph_error("*All* polygons undefined or out of range, thus no plot.");
  2737.     }
  2738.     sort_by_zmax();
  2739.  
  2740.     /* sort the polygons and draw them in one loop */
  2741.     Last = -1L;
  2742.     This = pfirst;
  2743.     /* Search first polygon, that can be drawn */
  2744.     while (Last == -1L && This >= 0L) {
  2745.     This = pfirst;
  2746.     if (obscured(plist + This)) {
  2747.         Last = This;
  2748.         continue;
  2749.     }
  2750.     if (plist[This].tested != is_tested && !in_front(Last, This))
  2751.         continue;
  2752.     draw(plist + This);
  2753.     Last = This;
  2754.     }
  2755.     /* Draw the other polygons */
  2756.     for (This = plist[Last].next; This >= 0L; This = plist[Last].next) {
  2757.     if (obscured(plist + This)) {
  2758.         Last = This;
  2759.         continue;
  2760.     }
  2761.     if (plist[This].tested != is_tested && !in_front(Last, This))
  2762.         continue;
  2763.     draw(plist + This);
  2764.     Last = This;
  2765.     }
  2766.  
  2767.     /* Free memory */
  2768.     if (plist) {
  2769.     for (This = 0L; This < pfree; This++)
  2770.         free(plist[This].vertex);
  2771.     free(plist);
  2772.     plist = 0;
  2773.     }
  2774.     if (vlist) {
  2775.     free(vlist);
  2776.     vlist = 0;
  2777.     }
  2778.     if (ymin_hl) {
  2779.     free(ymin_hl);
  2780.     ymin_hl = 0;
  2781.     }
  2782.     if (ymax_hl) {
  2783.     free(ymax_hl);
  2784.     ymax_hl = 0;
  2785.     }
  2786.     if (hl_buffer) {
  2787.     free(hl_buffer);
  2788.     hl_buffer = 0;
  2789.     }
  2790.     if (Cross_store) {
  2791.     free(Cross_store);
  2792.     Cross_store = 0;
  2793.     last_Cross_store = 0;
  2794.     }
  2795. }
  2796.